Changeset 3504 for trunk/LMDZ.PLUTO/libf
- Timestamp:
- Nov 8, 2024, 10:57:14 AM (6 weeks ago)
- Location:
- trunk/LMDZ.PLUTO/libf/phypluto
- Files:
-
- 10 edited
-
ave_stelspec.F90 (modified) (10 diffs)
-
callcorrk_pluto_mod.F90 (modified) (4 diffs)
-
iniorbit.F (modified) (2 diffs)
-
phyetat0_mod.F90 (modified) (12 diffs)
-
physiq_mod.F90 (modified) (1 diff)
-
setspi.F90 (modified) (14 diffs)
-
setspv.F90 (modified) (9 diffs)
-
suaer_corrk.F90 (modified) (1 diff)
-
sugas_corrk.F90 (modified) (35 diffs)
-
surfini.F90 (modified) (3 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.PLUTO/libf/phypluto/ave_stelspec.F90
r3184 r3504 2 2 3 3 !================================================================== 4 ! 4 ! 5 5 ! Purpose 6 6 ! ------- 7 7 ! Average the chosen high resolution stellar spectrum over the 8 8 ! visible bands in the model. 9 ! 9 ! 10 10 ! Authors 11 ! ------- 11 ! ------- 12 12 ! Robin Wordsworth (2010). 13 13 ! Generalized to very late spectral types (and Brown dwarfs) Jeremy Leconte (2012) … … 16 16 ! --------- 17 17 ! setspv.F 18 ! 18 ! 19 19 ! Calls 20 20 ! ----- 21 21 ! none 22 ! 22 ! 23 23 !================================================================== 24 24 … … 28 28 use callkeys_mod, only: stelbbody,stelTbb,startype 29 29 use ioipsl_getin_p_mod, only: getin_p 30 use mod_phys_lmdz_para, only : is_master 30 31 31 32 implicit none 32 33 33 34 real*8 STELLAR(L_NSPECTV) 34 ! integer startype 35 ! integer startype 35 36 36 37 integer Nfine … … 48 49 real lam_temp 49 50 double precision stel_temp 50 51 51 52 integer :: ios ! file opening/reading status 52 53 53 54 STELLAR(:)=0.0 54 55 55 print*,'enter ave_stellspec'56 if (is_master) print*,'enter ave_stellspec' 56 57 if(stelbbody)then 57 58 tstellar=stelTbb … … 65 66 call blackl(dble(lam_temp*1e-6),dble(tstellar),stel_temp) 66 67 STELLAR(band)=STELLAR(band)+stel_temp*dl 67 enddo 68 enddo 68 69 end do 69 70 STELLAR(1:L_NSPECTV)=STELLAR(1:L_NSPECTV)/sum(STELLAR(1:L_NSPECTV)) … … 99 100 Case(6) 100 101 file_id='/stellar_spectra/BD_Teff-1600K.txt' 101 tstellar=1600. 102 tstellar=1600. 102 103 file_id_lam='/stellar_spectra/lamBD.txt' 103 104 Nfine=5000 104 105 Case(7) 105 106 file_id='/stellar_spectra/BD_Teff-1000K.txt' 106 tstellar=1000. 107 tstellar=1000. 107 108 file_id_lam='/stellar_spectra/lamBD.txt' 108 109 Nfine=5000 109 110 Case(8) 110 111 file_id='/stellar_spectra/Flux_K5_Teff4700_logg4.5_Met-0.5_BTsettle.dat' 111 tstellar=4700. 112 tstellar=4700. 112 113 file_id_lam='/stellar_spectra/lambda_K5_Teff4700_logg4.5_Met-0.5_BTsettle.dat' 113 114 Nfine=3986 114 115 Case(9) 115 116 file_id='/stellar_spectra/Flux_TRAPPIST1.dat' 116 tstellar=2550. 117 tstellar=2550. 117 118 file_id_lam='/stellar_spectra/lambda_TRAPPIST1.dat' 118 119 Nfine=5000 119 120 Case(10) 120 121 file_id='/stellar_spectra/Flux_Proxima.dat' 121 tstellar=3050. 122 tstellar=3050. 122 123 file_id_lam='/stellar_spectra/lambda_Proxima.dat' 123 124 Nfine=5000 … … 139 140 stop 140 141 end if 141 142 142 143 ! Opening file 143 144 file_path = trim(datadir)//'/stellar_spectra/'//stelspec_file 144 145 print*, 'stellar flux : ', file_path 145 146 OPEN(UNIT=110,FILE=file_path,STATUS='old',iostat=ios) 146 147 147 148 if (ios /= 0) THEN 148 149 write(*,*)'Error: cannot open stelspec_file file ', trim(stelspec_file) … … 219 220 !$OMP END MASTER 220 221 !$OMP BARRIER 221 222 222 223 ! sum data by band 223 224 band=1 … … 236 237 STELLAR(band)=STELLAR(band)+stel_f(ifine)*dl 237 238 end do 238 239 239 240 240 241 STELLAR(1:L_NSPECTV)=STELLAR(1:L_NSPECTV)/sum(STELLAR(1:L_NSPECTV)) 241 242 !$OMP BARRIER … … 244 245 if (allocated(stel_f)) deallocate(stel_f) 245 246 !$OMP END MASTER 246 !$OMP BARRIER 247 !$OMP BARRIER 247 248 endif 248 249 -
trunk/LMDZ.PLUTO/libf/phypluto/callcorrk_pluto_mod.F90
r3501 r3504 32 32 use sfluxi_pluto_mod, only: sfluxi_pluto 33 33 use sfluxv_pluto_mod, only: sfluxv_pluto 34 use mod_phys_lmdz_para, only : is_master 34 35 35 36 … … 222 223 if(firstcall) then 223 224 224 print*, "callcorrk: Correlated-k data folder:",trim(datadir)225 if (is_master) print*, "callcorrk: Correlated-k data folder:",trim(datadir) 225 226 call getin("corrkdir",corrkdir) 226 227 print*, "corrkdir = ",corrkdir … … 230 231 banddir=trim(adjustl(corrkdir))//'/'//trim(adjustl(banddir)) 231 232 232 print*,'starting sugas'233 if (is_master) print*,'starting sugas' 233 234 call sugas_corrk ! set up gaseous absorption properties 234 print*,'starting setspi'235 if (is_master) print*,'starting setspi' 235 236 call setspi ! basic infrared properties 236 print*,'starting setspv'237 if (is_master) print*,'starting setspv' 237 238 call setspv ! basic visible properties 238 239 … … 259 260 call haze_reffrad_fix(ngrid,nlayer,zzlay, & 260 261 reffrad,nueffrad) 261 print*, 'haze_radproffix=T : fixed profile for haze rad'262 if (is_master) print*, 'haze_radproffix=T : fixed profile for haze rad' 262 263 else 263 print*,'reffrad haze:',reffrad(1,1,iaero_haze)264 print*,'nueff haze',nueffrad(1,1,iaero_haze)264 if (is_master) print*,'reffrad haze:',reffrad(1,1,iaero_haze) 265 if (is_master) print*,'nueff haze',nueffrad(1,1,iaero_haze) 265 266 endif 266 267 endif ! radiative haze -
trunk/LMDZ.PLUTO/libf/phypluto/iniorbit.F
r3184 r3504 1 1 SUBROUTINE iniorbit 2 2 $ (papoastr,pperiastr,pyear_day,pperi_day,pobliq) 3 3 4 4 USE planete_mod, only: apoastr, periastr, year_day, obliquit, 5 5 & peri_day, e_elips, p_elips, timeperi 6 6 use comcstfi_mod, only: pi 7 use mod_phys_lmdz_para, only : is_master 7 8 IMPLICIT NONE 8 9 … … 32 33 peri_day=pperi_day 33 34 34 PRINT*,'iniorbit: Periastron in AU ',periastr 35 PRINT*,'iniorbit: Apoastron in AU ',apoastr 36 PRINT*,'iniorbit: Obliquity in degrees :',obliquit 37 35 if (is_master) then 36 PRINT*,'iniorbit: Periastron in AU ',periastr 37 PRINT*,'iniorbit: Apoastron in AU ',apoastr 38 PRINT*,'iniorbit: Obliquity in degrees :',obliquit 39 endif 38 40 39 41 e_elips=(apoastr-periastr)/(periastr+apoastr) 40 42 p_elips=0.5*(periastr+apoastr)*(1-e_elips*e_elips) 41 43 42 print*,'iniorbit: e_elips',e_elips 43 print*,'iniorbit: p_elips',p_elips 44 44 if (is_master) then 45 print*,'iniorbit: e_elips',e_elips 46 print*,'iniorbit: p_elips',p_elips 47 endif 45 48 !----------------------------------------------------------------------- 46 49 ! compute polar angle and distance to the Sun: -
trunk/LMDZ.PLUTO/libf/phypluto/phyetat0_mod.F90
r3502 r3504 26 26 inquire_dimension, inquire_dimension_length 27 27 use callkeys_mod, only: surfalbedo,surfemis, callsoil 28 use mod_phys_lmdz_para, only : is_master 28 29 implicit none 29 30 … … 100 101 101 102 ! possibility to modify tab_cntrl in tabfi 102 write(*,*)103 write(*,*) 'TABFI in phyeta0: Lmodif=',Lmodif," tab0=",tab0103 if (is_master) write(*,*) 104 if (is_master) write(*,*) 'TABFI in phyeta0: Lmodif=',Lmodif," tab0=",tab0 104 105 call tabfi (ngrid,nid_start,Lmodif,tab0,day_ini,lmax,p_rad, & 105 106 p_omeg,p_g,p_cpp,p_mugaz,p_daysec,time) … … 119 120 phisfi(:)=0. 120 121 endif ! of if (startphy_file) 121 write(*,*) "phyetat0: surface geopotential <phisfi> range:", &122 if (is_master) write(*,*) "phyetat0: surface geopotential <phisfi> range:", & 122 123 minval(phisfi), maxval(phisfi) 123 124 … … 132 133 surfalbedo=0.5 133 134 call getin_p("surfalbedo",surfalbedo) 134 print*,"surfalbedo",surfalbedo135 if (is_master) print*,"surfalbedo",surfalbedo 135 136 albedodat(:)=surfalbedo 136 137 endif ! of if (startphy_file) 137 write(*,*) "phyetat0: Bare ground albedo <albedodat> range:", &138 if (is_master) write(*,*) "phyetat0: Bare ground albedo <albedodat> range:", & 138 139 minval(albedodat), maxval(albedodat) 139 140 … … 147 148 zmea(:)=0. 148 149 endif ! of if (startphy_file) 149 write(*,*) "phyetat0: <ZMEA> range:", &150 if (is_master) write(*,*) "phyetat0: <ZMEA> range:", & 150 151 minval(zmea), maxval(zmea) 151 152 … … 159 160 zstd(:)=0. 160 161 endif ! of if (startphy_file) 161 write(*,*) "phyetat0: <ZSTD> range:", &162 if (is_master) write(*,*) "phyetat0: <ZSTD> range:", & 162 163 minval(zstd), maxval(zstd) 163 164 … … 171 172 zsig(:)=0. 172 173 endif ! of if (startphy_file) 173 write(*,*) "phyetat0: <ZSIG> range:", &174 if (is_master) write(*,*) "phyetat0: <ZSIG> range:", & 174 175 minval(zsig), maxval(zsig) 175 176 … … 183 184 zgam(:)=0. 184 185 endif ! of if (startphy_file) 185 write(*,*) "phyetat0: <ZGAM> range:", &186 if (is_master) write(*,*) "phyetat0: <ZGAM> range:", & 186 187 minval(zgam), maxval(zgam) 187 188 … … 195 196 zthe(:)=0. 196 197 endif ! of if (startphy_file) 197 write(*,*) "phyetat0: <ZTHE> range:", &198 if (is_master) write(*,*) "phyetat0: <ZTHE> range:", & 198 199 minval(zthe), maxval(zthe) 199 200 … … 207 208 tsurf(:)=0. ! will be updated afterwards in physiq ! 208 209 endif ! of if (startphy_file) 209 write(*,*) "phyetat0: Surface temperature <tsurf> range:", &210 if (is_master) write(*,*) "phyetat0: Surface temperature <tsurf> range:", & 210 211 minval(tsurf), maxval(tsurf) 211 212 … … 220 221 surfemis=1.0 221 222 call getin_p("surfemis",surfemis) 222 print*,"surfemis",surfemis223 if (is_master) print*,"surfemis",surfemis 223 224 emis(:)=surfemis 224 225 endif ! of if (startphy_file) 225 write(*,*) "phyetat0: Surface emissivity <emis> range:", &226 if (is_master) write(*,*) "phyetat0: Surface emissivity <emis> range:", & 226 227 minval(emis), maxval(emis) 227 228 … … 237 238 q2(:,:)=0. 238 239 endif ! of if (startphy_file) 239 write(*,*) "phyetat0: PBL wind variance <q2> range:", &240 if (is_master) write(*,*) "phyetat0: PBL wind variance <q2> range:", & 240 241 minval(q2), maxval(q2) 241 242 -
trunk/LMDZ.PLUTO/libf/phypluto/physiq_mod.F90
r3501 r3504 726 726 ! at firstcall of callcorrk so we can output XspecIR, XspecVI 727 727 ! when using Dynamico 728 print*, "physiq_mod: Correlated-k data base folder:",trim(datadir)728 if (is_master) print*, "physiq_mod: Correlated-k data base folder:",trim(datadir) 729 729 call getin_p("corrkdir",corrkdir) 730 print*,"corrkdir = ", corrkdir730 if (is_master) print*,"corrkdir = ", corrkdir 731 731 write (tmp1, '(i4)') L_NSPECTI 732 732 write (tmp2, '(i4)') L_NSPECTV -
trunk/LMDZ.PLUTO/libf/phypluto/setspi.F90
r3184 r3504 2 2 3 3 !================================================================== 4 ! 4 ! 5 5 ! Purpose 6 6 ! ------- 7 7 ! Set up spectral intervals and Planck function in the longwave. 8 ! 8 ! 9 9 ! Authors 10 ! ------- 10 ! ------- 11 11 ! Adapted from setspi in the NASA Ames radiative code by 12 12 ! Robin Wordsworth (2009). 13 ! 13 ! 14 14 ! Called by 15 15 ! --------- 16 16 ! callcorrk.F 17 ! 17 ! 18 18 ! Calls 19 19 ! ----- 20 20 ! none 21 ! 21 ! 22 22 !================================================================== 23 23 … … 26 26 use datafile_mod, only: datadir 27 27 use comcstfi_mod, only: pi 28 use mod_phys_lmdz_para, only : is_master 28 29 29 30 implicit none … … 42 43 real*8 :: c1 = 3.741832D-16 ! W m^-2 43 44 real*8 :: c2 = 1.438786D-2 ! m K 44 45 45 46 real*8 :: lastband(2), plancksum 46 47 … … 75 76 write(temp1,'(i2.2)') L_NSPECTI 76 77 !file_id='/corrk_data/' // corrkdir(1:LEN_TRIM(corrkdir)) // '/narrowbands_IR.in' 77 file_id='/corrk_data/'//trim(adjustl(banddir))//'/narrowbands_IR.in' 78 file_id='/corrk_data/'//trim(adjustl(banddir))//'/narrowbands_IR.in' 78 79 file_path=TRIM(datadir)//TRIM(file_id) 79 80 … … 81 82 inquire(FILE=file_path,EXIST=file_ok) 82 83 if(.not.file_ok) then 83 write(*,*)'The file ',TRIM(file_path) 84 write(*,*)'was not found by setspi.F90, exiting.' 85 write(*,*)'Check that your path to datagcm:',trim(datadir) 86 write(*,*)' is correct. You can change it in callphys.def with:' 87 write(*,*)' datadir = /absolute/path/to/datagcm' 88 write(*,*)'Also check that the corrkdir you chose in callphys.def exists.' 89 call abort 90 endif 91 92 !$OMP MASTER 84 if (is_master) then 85 write(*,*)'The file ',TRIM(file_path) 86 write(*,*)'was not found by setspi.F90, exiting.' 87 write(*,*)'Check that your path to datagcm:',trim(datadir) 88 write(*,*)' is correct. You can change it in callphys.def with:' 89 write(*,*)' datadir = /absolute/path/to/datagcm' 90 write(*,*)'Also check that the corrkdir you chose in callphys.def exists.' 91 endif 92 call abort_physic("setspi", "File not found by setspi", 1) 93 endif 94 95 !$OMP MASTER 93 96 nb=0 94 97 ierr=0 95 ! check that the file contains the right number of bands 98 ! check that the file contains the right number of bands 96 99 open(131,file=file_path,form='formatted') 97 100 read(131,*,iostat=ierr) file_entries … … 103 106 close(131) 104 107 105 write(*,*) 'setspi: L_NSPECTI = ',L_NSPECTI, 'in the model ' 106 write(*,*) ' there are ',nb, 'entries in ',TRIM(file_path) 108 if (is_master) then 109 write(*,*) 'setspi: L_NSPECTI = ',L_NSPECTI, 'in the model ' 110 write(*,*) ' there are ',nb, 'entries in ',TRIM(file_path) 111 endif 107 112 if(nb.ne.L_NSPECTI) then 108 write(*,*) 'MISMATCH !! I stop here' 109 call abort 113 call abort_physic("setspi",'MISMATCH !! I stop here',1) 110 114 endif 111 115 112 116 ! load and display the data 113 117 open(111,file=file_path,form='formatted') 114 read(111,*) 118 read(111,*) 115 119 do M=1,L_NSPECTI-1 116 120 read(111,*) BWNI(M) … … 123 127 !$OMP BARRIER 124 128 125 print*,'' 126 print*,'setspi: IR band limits:' 127 do M=1,L_NSPECTI+1 128 print*,m,'-->',BWNI(M),' cm^-1' 129 end do 130 131 ! Set up mean wavenumbers and wavenumber deltas. Units of 129 if (is_master)then 130 print*,'' 131 print*,'setspi: IR band limits:' 132 do M=1,L_NSPECTI+1 133 print*,m,'-->',BWNI(M),' cm^-1' 134 end do 135 endif 136 137 ! Set up mean wavenumbers and wavenumber deltas. Units of 132 138 ! wavenumbers is cm^(-1); units of wavelengths is microns. 133 139 … … 136 142 DWNI(M) = BWNI(M+1)-BWNI(M) 137 143 WAVEI(M) = 1.0D+4/WNOI(M) 138 BLAMI(M) = 0.01D0/BWNI(M) 144 BLAMI(M) = 0.01D0/BWNI(M) 139 145 end do 140 146 BLAMI(M) = 0.01D0/BWNI(M) … … 147 153 ! original planck.f; W m^-2 wavenumber^-1, where wavenumber is in CM^-1. 148 154 149 print*,'' 150 print*,'setspi: Current Planck integration range:' 151 print*,'T = ',dble(NTstart)/NTfac, ' to ',dble(NTstop)/NTfac,' K.' 155 if (is_master)then 156 print*,'' 157 print*,'setspi: Current Planck integration range:' 158 print*,'T = ',dble(NTstart)/NTfac, ' to ',dble(NTstop)/NTfac,' K.' 159 endif 152 160 153 161 IF(.NOT.ALLOCATED(planckir)) ALLOCATE(planckir(L_NSPECTI,NTstop-NTstart+1)) … … 159 167 bma = (b-a)/2.0D0 160 168 ! if (nw .eq. 25) then !LT debug 161 ! print*, "a = ",a 162 ! print*, "b= ",b 163 ! print*,"bpa = ",bpa 169 ! print*, "a = ",a 170 ! print*, "b= ",b 171 ! print*,"bpa = ",bpa 164 172 ! print*, "bma = ",bma 165 173 ! endif … … 170 178 y = bma*x(mm)+bpa 171 179 !to avoid floating overflow when T is low and optical wavelength 172 if ((c2/(y*T)) .lt. 700.0D0) then 180 if ((c2/(y*T)) .lt. 700.0D0) then 173 181 ans = ans + w(mm)*c1/(y**5*(exp(c2/(y*T))-1.0D0)) 174 else 182 else 175 183 ans = ans +0.0D0 176 184 endif … … 179 187 end do 180 188 end do 181 189 182 190 ! force planck=sigma*eps*T^4 for each temperature in array 183 191 if(forceEC)then 184 print*,'setspi: Force F=sigma*eps*T^4 for all values of T!'192 if (is_master) print*,'setspi: Force F=sigma*eps*T^4 for all values of T!' 185 193 do nt=NTstart,NTstop 186 194 plancksum=0.0D0 187 195 T=dble(NT)/NTfac 188 196 189 197 do NW=1,L_NSPECTI 190 198 plancksum=plancksum+ & … … 207 215 plancksum=plancksum+planckir(NW,nt-NTstart+1)*DWNI(NW)*pi 208 216 end do 209 print*,'setspi: At lower limit:' 210 print*,'in model sig*T^4 = ',plancksum,' W m^-2' 211 print*,'actual sig*T^4 = ',sigma*(dble(nt)/NTfac)**4,' W m^-2' 212 217 if (is_master) then 218 print*,'setspi: At lower limit:' 219 print*,'in model sig*T^4 = ',plancksum,' W m^-2' 220 print*,'actual sig*T^4 = ',sigma*(dble(nt)/NTfac)**4,' W m^-2' 221 endif 213 222 ! check energy conservation at upper temperature boundary 214 223 plancksum=0.0D0 … … 217 226 plancksum=plancksum+planckir(NW,nt-NTstart+1)*DWNI(NW)*pi 218 227 end do 219 print*,'setspi: At upper limit:' 220 print*,'in model sig*T^4 = ',plancksum,' W m^-2' 221 print*,'actual sig*T^4 = ',sigma*(dble(nt)/NTfac)**4,' W m^-2' 222 print*,'' 228 if (is_master) then 229 print*,'setspi: At upper limit:' 230 print*,'in model sig*T^4 = ',plancksum,' W m^-2' 231 print*,'actual sig*T^4 = ',sigma*(dble(nt)/NTfac)**4,' W m^-2' 232 print*,'' 233 endif 223 234 endif 224 235 -
trunk/LMDZ.PLUTO/libf/phypluto/setspv.F90
r3184 r3504 2 2 3 3 !================================================================== 4 ! 4 ! 5 5 ! Purpose 6 6 ! ------- 7 ! Set up spectral intervals, stellar spectrum and Rayleigh 8 ! opacity in the shortwave. 9 ! 7 ! Set up spectral intervals, stellar spectrum and Rayleigh 8 ! opacity in the shortwave. 9 ! 10 10 ! Authors 11 ! ------- 11 ! ------- 12 12 ! Adapted from setspv in the NASA Ames radiative code by 13 13 ! Robin Wordsworth (2009). … … 16 16 ! --------- 17 17 ! callcorrk.F 18 ! 18 ! 19 19 ! Calls 20 20 ! ----- 21 21 ! ave_stelspec.F 22 ! 22 ! 23 23 !================================================================== 24 24 … … 28 28 use datafile_mod, only: datadir 29 29 use callkeys_mod, only: Fat1AU,rayleigh 30 use mod_phys_lmdz_para, only : is_master 30 31 31 32 implicit none … … 53 54 54 55 write(temp1,'(i2.2)') L_NSPECTV 55 file_id='/corrk_data/'//trim(adjustl(banddir))//'/narrowbands_VI.in' 56 file_id='/corrk_data/'//trim(adjustl(banddir))//'/narrowbands_VI.in' 56 57 file_path=TRIM(datadir)//TRIM(file_id) 57 58 … … 65 66 write(*,*)' datadir = /absolute/path/to/datagcm' 66 67 write(*,*)'Also check that the corrkdir you chose in callphys.def exists.' 67 call abort 68 call abort_physic("setspv", "File not found by setspi", 1) 68 69 endif 69 70 !$OMP MASTER 70 71 !$OMP MASTER 71 72 nb=0 72 73 ierr=0 73 ! check that the file contains the right number of bands 74 ! check that the file contains the right number of bands 74 75 open(131,file=file_path,form='formatted') 75 76 read(131,*,iostat=ierr) file_entries … … 80 81 close(131) 81 82 82 write(*,*) 'setspv: L_NSPECTV = ',L_NSPECTV, 'in the model ' 83 write(*,*) ' there are ',nb, 'entries in ',TRIM(file_path) 83 if (is_master) then 84 write(*,*) 'setspv: L_NSPECTV = ',L_NSPECTV, 'in the model ' 85 write(*,*) ' there are ',nb, 'entries in ',TRIM(file_path) 86 endif 84 87 if(nb.ne.L_NSPECTV) then 85 write(*,*) 'MISMATCH !! I stop here' 86 call abort 88 call abort_physic("setspv",'MISMATCH !! I stop here',1) 87 89 endif 88 90 89 91 ! load and display the data 90 92 open(111,file=file_path,form='formatted') 91 read(111,*) 93 read(111,*) 92 94 do M=1,L_NSPECTV-1 93 95 read(111,*) BWNV(M) … … 100 102 !$OMP BARRIER 101 103 102 print*,'setspv: VI band limits:' 103 do M=1,L_NSPECTV+1 104 print*,m,'-->',BWNV(M),' cm^-1' 105 end do 106 print*,' ' 104 if (is_master) then 105 print*,'setspv: VI band limits:' 106 do M=1,L_NSPECTV+1 107 print*,m,'-->',BWNV(M),' cm^-1' 108 end do 109 print*,' ' 110 endif 107 111 108 ! Set up mean wavenumbers and wavenumber deltas. Units of 112 ! Set up mean wavenumbers and wavenumber deltas. Units of 109 113 ! wavenumbers is cm^(-1); units of wavelengths is microns. 110 114 … … 121 125 ! Set up stellar spectrum 122 126 123 write(*,*)'setspv: Interpolating stellar spectrum from the hires data...'127 if (is_master) write(*,*)'setspv: Interpolating stellar spectrum from the hires data...' 124 128 call ave_stelspec(STELLAR) 125 129 126 ! Sum the stellar flux, and write out the result. 127 sum = 0.0 130 ! Sum the stellar flux, and write out the result. 131 sum = 0.0 128 132 do N=1,L_NSPECTV 129 133 STELLARF(N) = STELLAR(N) * Fat1AU 130 134 sum = sum+STELLARF(N) 131 135 end do 132 write(6,'("setspv: Stellar flux at 1 AU = ",f9.2," W m-2")') sum 133 print*,' ' 136 if (is_master)then 137 write(6,'("setspv: Stellar flux at 1 AU = ",f9.2," W m-2")') sum 138 print*,' ' 139 endif 134 140 135 141 … … 142 148 call calc_rayleigh 143 149 else 144 print*,'setspv: No Rayleigh scattering, check for NaN in output!'150 if (is_master) print*,'setspv: No Rayleigh scattering, check for NaN in output!' 145 151 do N=1,L_NSPECTV 146 152 TAURAY(N) = 1E-16 -
trunk/LMDZ.PLUTO/libf/phypluto/suaer_corrk.F90
r3501 r3504 142 142 143 143 ! Only one table of optical properties : 144 write(*,*)'Suaer haze optical properties, using: ', &144 if (is_master) write(*,*)'Suaer haze optical properties, using: ', & 145 145 TRIM(hazeprop_file) 146 146 -
trunk/LMDZ.PLUTO/libf/phypluto/sugas_corrk.F90
r3184 r3504 38 38 use recombin_corrk_mod, only: su_recombin, & 39 39 corrk_recombin, use_premix, nrecomb_tot, rcb2tot_idx 40 use mod_phys_lmdz_para, only : is_master 40 41 implicit none 41 42 … … 65 66 66 67 if (.not. moderntracdef) use_premix=.true. ! Added by JVO for compatibility with 'old' traceur.def 67 68 68 69 !$OMP MASTER 69 70 if (use_premix) then ! use_premix flag added by JVO, thus if pure recombining then premix is skipped … … 98 99 'match that in gases.def (',ngasmx,'), exiting.' 99 100 call abort 100 endif 101 endif 101 102 endif 102 103 … … 111 112 '], radiative code cannot handle this.' 112 113 call abort 113 endif 114 endif 114 115 115 116 ! dynamically allocate gastype and read from Q.dat … … 145 146 ! Check that gastype and gnom match 146 147 do igas=1,ngas 147 print*,'Gas ',igas,' is ',trim(gnom(igas))148 if (is_master) print*,'Gas ',igas,' is ',trim(gnom(igas)) 148 149 if (trim(gnom(igas)).ne.trim(gastype(igas))) then 149 150 print*,'Name of a gas in radiative transfer data (',trim(gastype(igas)),') does not ', & 150 151 'match that in gases.def (',trim(gnom(igas)),'), exiting. You should compare ', & 151 'gases.def with Q.dat in your radiative transfer directory.' 152 'gases.def with Q.dat in your radiative transfer directory.' 152 153 call abort 153 154 endif 154 155 enddo 155 print*,'Confirmed gas match in radiative transfer and gases.def!'156 if (is_master) print*,'Confirmed gas match in radiative transfer and gases.def!' 156 157 endif 157 158 158 159 ! display the values 159 print*,'Variable gas volume mixing ratios:' 160 do n=1,L_REFVAR 161 !print*,n,'.',wrefvar(n),' kg/kg' ! pay attention! 162 print*,n,'.',wrefvar(n),' mol/mol' 163 end do 164 print*,'' 165 160 if (is_master)then 161 print*,'Variable gas volume mixing ratios:' 162 do n=1,L_REFVAR 163 !print*,n,'.',wrefvar(n),' kg/kg' ! pay attention! 164 print*,n,'.',wrefvar(n),' mol/mol' 165 end do 166 print*,'' 167 endif 168 166 169 else 167 170 L_REFVAR = 1 … … 185 188 call abort 186 189 endif 187 190 188 191 ! check the array size is correct, load the coefficients 189 192 open(111,file=TRIM(file_path),form='formatted') … … 192 195 read(111,*) gweight 193 196 close(111) 194 197 195 198 ! display the values 196 print*,'Correlated-k g-space grid:' 197 do n=1,L_NGAUSS 198 print*,n,'.',gweight(n) 199 end do 200 print*,'' 199 if (is_master)then 200 print*,'Correlated-k g-space grid:' 201 do n=1,L_NGAUSS 202 print*,n,'.',gweight(n) 203 end do 204 print*,'' 205 endif 201 206 202 207 !======================================================================= … … 221 226 call abort 222 227 endif 223 228 224 229 ! get array size, load the coefficients 225 230 open(111,file=TRIM(file_path),form='formatted') … … 232 237 233 238 ! display the values 234 print*,'Correlated-k pressure grid (mBar):' 235 do n=1,L_NPREF 236 print*,n,'. 1 x 10^',pgasref(n),' mBar' 237 end do 238 print*,'' 239 if (is_master)then 240 print*,'Correlated-k pressure grid (mBar):' 241 do n=1,L_NPREF 242 print*,n,'. 1 x 10^',pgasref(n),' mBar' 243 end do 244 print*,'' 245 endif 239 246 240 247 ! save the min / max matrix values … … 276 283 277 284 ! display the values 278 print*,'Correlated-k temperature grid:' 279 do n=1,L_NTREF 280 print*,n,'.',tgasref(n),' K' 281 end do 285 if (is_master)then 286 print*,'Correlated-k temperature grid:' 287 do n=1,L_NTREF 288 print*,n,'.',tgasref(n),' K' 289 end do 290 endif 282 291 283 292 ! save the min / max matrix values … … 299 308 !----------------------------------------------------------------------- 300 309 ! allocate the multidimensional arrays in radcommon_h 301 if(corrk_recombin) then 302 IF(.NOT. ALLOCATED(gasi)) ALLOCATE(gasi(L_NTREF,L_PINT,L_REFVAR+nrecomb_tot,L_NSPECTI,L_NGAUSS)) 303 IF(.NOT. ALLOCATED(gasv)) ALLOCATE(gasv(L_NTREF,L_PINT,L_REFVAR+nrecomb_tot,L_NSPECTV,L_NGAUSS)) 304 ! display the values 305 print*,'' 306 print*,'Correlated-k matrix size:' 307 print*,'[',L_NTREF,',',L_NPREF,',',L_REFVAR+nrecomb_tot,' (',L_REFVAR,'+',nrecomb_tot,'),',L_NGAUSS,']' 310 if(corrk_recombin) then 311 IF(.NOT. ALLOCATED(gasi)) ALLOCATE(gasi(L_NTREF,L_PINT,L_REFVAR+nrecomb_tot,L_NSPECTI,L_NGAUSS)) 312 IF(.NOT. ALLOCATED(gasv)) ALLOCATE(gasv(L_NTREF,L_PINT,L_REFVAR+nrecomb_tot,L_NSPECTV,L_NGAUSS)) 313 ! display the values 314 if (is_master)then 315 print*,'' 316 print*,'Correlated-k matrix size:' 317 print*,'[',L_NTREF,',',L_NPREF,',',L_REFVAR+nrecomb_tot,' (',L_REFVAR,'+',nrecomb_tot,'),',L_NGAUSS,']' 318 endif 308 319 else 309 IF(.NOT. ALLOCATED(gasi)) ALLOCATE(gasi(L_NTREF,L_PINT,L_REFVAR,L_NSPECTI,L_NGAUSS)) 310 IF(.NOT. ALLOCATED(gasv)) ALLOCATE(gasv(L_NTREF,L_PINT,L_REFVAR,L_NSPECTV,L_NGAUSS)) 311 ! display the values 312 print*,'' 313 print*,'Correlated-k matrix size:' 314 print*,'[',L_NTREF,',',L_NPREF,',',L_REFVAR,',',L_NGAUSS,']' 320 IF(.NOT. ALLOCATED(gasi)) ALLOCATE(gasi(L_NTREF,L_PINT,L_REFVAR,L_NSPECTI,L_NGAUSS)) 321 IF(.NOT. ALLOCATED(gasv)) ALLOCATE(gasv(L_NTREF,L_PINT,L_REFVAR,L_NSPECTV,L_NGAUSS)) 322 ! display the values 323 if (is_master)then 324 print*,'' 325 print*,'Correlated-k matrix size:' 326 print*,'[',L_NTREF,',',L_NPREF,',',L_REFVAR,',',L_NGAUSS,']' 327 endif 315 328 endif 316 329 … … 321 334 ! wavelength used to separate IR from VI in graybody. We will need that anyway 322 335 IR_VI_wnlimit=3000. 323 write(*,*)"graybody: Visible / Infrared separation set at",10000./IR_VI_wnlimit,"um"324 336 if (is_master)write(*,*)"graybody: Visible / Infrared separation set at",10000./IR_VI_wnlimit,"um" 337 325 338 nVI_limit=0 326 339 Do nw=1,L_NSPECTV … … 344 357 call getin_p("kappa_VI",kappa_VI) 345 358 write(*,*)" kappa_VI = ",kappa_VI 346 kappa_VI=kappa_VI*1.e4* mugaz * 1.672621e-27 ! conversion from m^2/kg to cm^2/molecule 347 359 kappa_VI=kappa_VI*1.e4* mugaz * 1.672621e-27 ! conversion from m^2/kg to cm^2/molecule 360 348 361 ! constant absorption coefficient in IR 349 362 write(*,*)"graybody: constant absorption coefficient in InfraRed:" 350 363 kappa_IR=-100000. 351 364 call getin_p("kappa_IR",kappa_IR) 352 write(*,*)" kappa_IR = ",kappa_IR 353 kappa_IR=kappa_IR*1.e4* mugaz * 1.672621e-27 ! conversion from m^2/kg to cm^2/molecule 365 write(*,*)" kappa_IR = ",kappa_IR 366 kappa_IR=kappa_IR*1.e4* mugaz * 1.672621e-27 ! conversion from m^2/kg to cm^2/molecule 354 367 355 368 write(*,*)"graybody: Visible / Infrared separation set at band: IR=",nIR_limit,", VI=",nVI_limit 356 369 357 370 Else 358 kappa_VI=1.e-30 359 kappa_IR=1.e-30 371 kappa_VI=1.e-30 372 kappa_IR=1.e-30 360 373 End if 361 374 362 !$OMP MASTER 375 !$OMP MASTER 363 376 364 377 ! VISIBLE 365 378 if (callgasvis) then 366 379 367 380 ! Looking for premixed corrk files corrk_gcm_VI.dat if needed 368 381 if (use_premix) then … … 373 386 print*,'Visible corrk gaseous absorption is set to zero if graybody=F' 374 387 else 375 file_id='/corrk_data/'//trim(adjustl(banddir))//'/corrk_gcm_VI.dat' 388 file_id='/corrk_data/'//trim(adjustl(banddir))//'/corrk_gcm_VI.dat' 376 389 file_path=TRIM(datadir)//TRIM(file_id) 377 390 378 391 ! check that the file exists 379 392 inquire(FILE=file_path,EXIST=file_ok) … … 384 397 call abort 385 398 endif 386 399 387 400 open(111,file=TRIM(file_path),form='formatted') 388 401 read(111,*) gasv8(:,:,:L_REFVAR,:,:) … … 392 405 gasv8(:,:,1,:,:)=0.0 ! dummy but needs to be initialized 393 406 endif ! use_premix 394 407 395 408 ! Looking for others files corrk_gcm_VI_XXX.dat if needed 396 409 if (corrk_recombin) then 397 410 do igas=1,nrecomb_tot 398 411 399 file_id='/corrk_data/'//trim(adjustl(banddir))//'/corrk_gcm_VI_'//trim(adjustl(noms(rcb2tot_idx(igas))))//'.dat' 412 file_id='/corrk_data/'//trim(adjustl(banddir))//'/corrk_gcm_VI_'//trim(adjustl(noms(rcb2tot_idx(igas))))//'.dat' 400 413 file_path=TRIM(datadir)//TRIM(file_id) 401 414 402 415 ! check that the file exists 403 416 inquire(FILE=file_path,EXIST=file_ok) … … 408 421 call abort 409 422 endif 410 423 411 424 open(111,file=TRIM(file_path),form='formatted') 412 425 read(111,*) gasv8(:,:,L_REFVAR+igas,:,:) 413 426 close(111) 414 enddo 427 enddo 415 428 endif ! corrk_recombin 416 429 … … 423 436 gasv8(:,:,:,nVI_limit+1:L_NSPECTV,:)= gasv8(:,:,:,nVI_limit+1:L_NSPECTV,:)+kappa_VI 424 437 end if 425 426 else 438 439 else 427 440 print*,'Visible corrk gaseous absorption is set to zero.' 428 441 gasv8(:,:,:,:,:)=0.0 429 442 endif ! callgasvis 430 443 431 444 !$OMP END MASTER 432 445 !$OMP BARRIER 433 446 434 447 ! INFRA-RED 435 448 436 449 ! Looking for premixed corrk files corrk_gcm_IR.dat if needed 437 450 if (use_premix) then 438 451 if ((corrkdir(1:4).eq.'null'))then !.or.(TRIM(corrkdir).eq.'null_LowTeffStar')) then 439 452 print*,'Infrared corrk gaseous absorption is set to zero if graybody=F' 440 !$OMP MASTER 453 !$OMP MASTER 441 454 gasi8(:,:,:,:,:)=0.0 442 455 !$OMP END MASTER 443 456 !$OMP BARRIER 444 else 445 file_id='/corrk_data/'//trim(adjustl(banddir))//'/corrk_gcm_IR.dat' 457 else 458 file_id='/corrk_data/'//trim(adjustl(banddir))//'/corrk_gcm_IR.dat' 446 459 file_path=TRIM(datadir)//TRIM(file_id) 447 460 448 461 ! check that the file exists 449 462 inquire(FILE=file_path,EXIST=file_ok) … … 454 467 call abort 455 468 endif 456 457 !$OMP MASTER 469 470 !$OMP MASTER 458 471 open(111,file=TRIM(file_path),form='formatted') 459 472 read(111,*) gasi8(:,:,:L_REFVAR,:,:) … … 461 474 !$OMP END MASTER 462 475 !$OMP BARRIER 463 476 464 477 ! 'fzero' is a currently unused feature that allows optimisation 465 478 ! of the radiative transfer by neglecting bands where absorption 466 ! is close to zero. As it could be useful in the future, this 479 ! is close to zero. As it could be useful in the future, this 467 480 ! section of the code has been kept commented and not erased. 468 481 ! RW 7/3/12. … … 505 518 gasi8(:,:,1,:,:)=0.0 ! dummy but needs to be initialized 506 519 endif ! use_premix 507 520 508 521 ! Looking for others files corrk_gcm_IR_XXX.dat if needed 509 522 if (corrk_recombin) then 510 523 do igas=1,nrecomb_tot 511 524 512 file_id='/corrk_data/'//trim(adjustl(banddir))//'/corrk_gcm_IR_'//trim(adjustl(noms(rcb2tot_idx(igas))))//'.dat' 525 file_id='/corrk_data/'//trim(adjustl(banddir))//'/corrk_gcm_IR_'//trim(adjustl(noms(rcb2tot_idx(igas))))//'.dat' 513 526 file_path=TRIM(datadir)//TRIM(file_id) 514 527 515 528 ! check that the file exists 516 529 inquire(FILE=file_path,EXIST=file_ok) … … 521 534 call abort 522 535 endif 523 !$OMP MASTER 536 !$OMP MASTER 524 537 open(111,file=TRIM(file_path),form='formatted') 525 538 read(111,*) gasi8(:,:,L_REFVAR+igas,:,:) … … 527 540 !$OMP END MASTER 528 541 !$OMP BARRIER 529 enddo 542 enddo 530 543 endif ! corrk_recombin 531 544 532 !$OMP MASTER 545 !$OMP MASTER 533 546 if(nIR_limit.eq.0) then 534 547 gasi8(:,:,:,:,:)= gasi8(:,:,:,:,:)+kappa_VI … … 564 577 end if 565 578 end do 566 579 567 580 end do 568 581 end do … … 581 594 ! First, the initial interval 582 595 583 n = 1 596 n = 1 584 597 do m=1,5 585 598 x = pfgasref(m) … … 594 607 call lagrange(x,xi,yi,ans) 595 608 gasi(nt,m,nh,nw,ng) = 10.0**ans 596 end do 597 609 end do 610 598 611 do n=2,L_NPREF-2 599 612 do m=1,5 … … 610 623 call lagrange(x,xi,yi,ans) 611 624 gasi(nt,i,nh,nw,ng) = 10.0**ans 612 end do 625 end do 613 626 end do 614 627 615 628 ! Now, get the last interval 616 629 617 n = L_NPREF-1 630 n = L_NPREF-1 618 631 do m=1,5 619 632 i = (n-1)*5+m … … 629 642 call lagrange(x,xi,yi,ans) 630 643 gasi(nt,i,nh,nw,ng) = 10.0**ans 631 end do 644 end do 632 645 633 646 ! Fill the last pressure point … … 650 663 ! First, the initial interval 651 664 652 n = 1 665 n = 1 653 666 do m=1,5 654 667 x = pfgasref(m) … … 663 676 call lagrange(x,xi,yi,ans) 664 677 gasv(nt,m,nh,nw,ng) = 10.0**ans 665 end do 666 678 end do 679 667 680 do n=2,L_NPREF-2 668 681 do m=1,5 … … 679 692 call lagrange(x,xi,yi,ans) 680 693 gasv(nt,i,nh,nw,ng) = 10.0**ans 681 end do 694 end do 682 695 end do 683 696 … … 698 711 call lagrange(x,xi,yi,ans) 699 712 gasv(nt,i,nh,nw,ng) = 10.0**ans 700 end do 713 end do 701 714 702 715 ! Fill the last pressure point … … 704 717 gasv(nt,L_PINT,nh,nw,ng) = & 705 718 10.0**gasv8(nt,L_NPREF,nh,nw,ng) 706 719 707 720 end do 708 721 end do … … 743 756 ! endif 744 757 ! enddo 745 758 746 759 ! elseif (igas .eq. igas_H2O) then !AF24: removed 747 760 … … 750 763 endif 751 764 752 print*,'----------------------------------------------------' 753 print*,'And that`s all we have. It`s possible that other' 754 print*,'continuum absorption may be present, but if it is we' 755 print*,'don`t yet have data for it...' 756 print*,'' 765 if (is_master)then 766 print*,'----------------------------------------------------' 767 print*,'And that`s all we have. It`s possible that other' 768 print*,'continuum absorption may be present, but if it is we' 769 print*,'don`t yet have data for it...' 770 print*,'' 771 endif 757 772 758 773 ! Deallocate local arrays -
trunk/LMDZ.PLUTO/libf/phypluto/surfini.F90
r3503 r3504 7 7 use radinc_h, only : L_NSPECTV 8 8 use callkeys_mod, only : albedosnow, albedon2ice 9 use mod_phys_lmdz_para, only : is_master 9 10 10 11 IMPLICIT NONE … … 51 52 call planetwide_minval(albedo_bareground,min_albedo) 52 53 call planetwide_maxval(albedo_bareground,max_albedo) 53 write(*,*) 'surfini: minimum bare ground albedo',min_albedo54 write(*,*) 'surfini: maximum bare ground albedo',max_albedo54 if (is_master) write(*,*) 'surfini: minimum bare ground albedo',min_albedo 55 if (is_master) write(*,*) 'surfini: maximum bare ground albedo',max_albedo 55 56 56 57 … … 65 66 ENDDO 66 67 else 67 write(*,*) "surfini: No N2 ice tracer on surface ..."68 write(*,*) " and therefore no albedo change."68 if (is_master) write(*,*) "surfini: No N2 ice tracer on surface ..." 69 if (is_master) write(*,*) " and therefore no albedo change." 69 70 endif 70 71 call planetwide_minval(albedo,min_albedo) 71 72 call planetwide_maxval(albedo,max_albedo) 72 write(*,*) 'surfini: minimum corrected initial albedo',min_albedo73 write(*,*) 'surfini: maximum corrected initial albedo',max_albedo73 if (is_master) write(*,*) 'surfini: minimum corrected initial albedo',min_albedo 74 if (is_master) write(*,*) 'surfini: maximum corrected initial albedo',max_albedo 74 75 75 76
Note: See TracChangeset
for help on using the changeset viewer.
