! module init_module use params_module,only: nlon,nlat,nlevp1,nlonp4,nmlonp1,nmlat, | nmlev,glon1,dlon,glat1,dlat,nlev,plev1,dlev implicit none ! ! Initialize and store model non-input variables and constants: ! real :: glon(nlon),glat(nlat),plev(nlevp1) real :: gmlon(nmlonp1),gmlat(nmlat),pmlev(nmlev) ! ! istep: the current time step index. istep is initialized ! to 0, then incremented before each time step (see advance.f). ! integer :: istep ! time step index ! ! iter (iteration number): the number of timesteps (at the current ! step length) from model time 0,0,0 to current model time (day,hour,min). ! integer :: iter ! ! iyear and iday represent the current model calendar date ! (iyear is 4-digit). Uthr is current decimal hours. These ! are updated per timestep in advnce. ! integer :: iyear ! Current model calendar 4-digit year integer :: iday ! Current model calendar day real :: uthr ! Current ut (decimal hours) integer :: | start_mtime(3) ! starting model time (day,hr,min) ! ! getgpi flag must go here rather than in gpi_mod to avoid circular ! module dependency with nchist_mod and init_mod. ! integer :: igetgpi ! 0/1 flag to get GPI data (see gpi_ncfile ! in input_mod.f. Set in init_mod.f) integer :: igetgswmdi ! 0/1 flag to get GSWM data diurnal tides !(see gswm_di_ncfilein input_mod.f. Set in init_mod.f) integer :: igetgswmsdi! 0/1 flag to get GSWM data semidiurnal tides !(see gswm_sdi_ncfilein input_mod.f. Set in init_mod.f) integer :: igetgswmnmidi ! 0/1 flag to get GSWM data nonmigrating diurnal tides !(see gswm_nmidi_ncfilein input_mod.f. Set in init_mod.f) integer :: igetgswmnmisdi! 0/1 flag to get GSWM data nonmigrating semidiurnal tides !(see gswm_nmisdi_ncfilein input_mod.f. Set in init_mod.f) ! character(len=16) :: | host, ! host machine | system, ! operating system of host (from pre-proc macros) | logname ! user login name character(len=8) :: | rundate, ! current local date of run | runtime ! current local time of run ! ! Secs is updated in advnce, and is referenced in opflux, settei, ! sun, and chapmn. real :: | secs, ! current ut in seconds | sfeps, ! solar output change due to orbital eccentricity ! | alfalp,efluxlp ! low-energy protons in aurora | sundec, ! sun's declination (updated in advnce) | sin_sundec, ! sin(sundec) (updated in advnceday) | cos_sundec ! cos(sundec) (updated in advnceday) ! ! Day/night index is set by chapman.F: integer :: idn(nlonp4) ! day/night index ! contains !----------------------------------------------------------------------- subroutine init ! ! Initialize (this is called by tgcm.F after input): ! Some init also takes place in inp_model (input_mod.f) ! use input_module,only: start,step,secflds,secfmag,secfmagphr, | start_year,start_day,calendar_advance,gpi_ncfile,gswm_di_ncfile, | gswm_sdi_ncfile,gswm_nmidi_ncfile,gswm_nmisdi_ncfile, | mxhist_prim,mxhist_sech,output,secout,mkhvols,nmc,source_start use hist_module,only: hist_init,isechist,nfsech_geo,nfsech_mag, | nfsech_geo2d,nfsech_mag2d,nfsech_magphr,nstep,nhist_total, | nsech_total,nsource,nseries_prim,nseries_sech,nfiles_prim, | nfiles_sech use fields_module,only: set_fsech,init_fields, | fsech,fsechmag,fsech2d,fsechmag2d,fsechmagphr2d use cons_module,only: pi,init_cons,rtd,ylonm,ylatm use filter_module,only: trigs,ifax use mpi_module,only: lon0,lon1,lat0,lat1,mytid ! ! External: integer,external :: mtime_to_nstep ! ! Local: real :: theta0 integer :: i,iprintf ! #if (MPI == 0) lon0 = 1 lon1 = nlonp4 lat0 = 1 lat1 = nlat #endif ! ! Initialize derived constants (init_cons is in cons_module): call init_cons ! ! Mag grid coordinates ylonm,ylatm were set in init_consdyn (cons.F) ! (convert to degrees): do i=1,nmlonp1 gmlon(i) = ylonm(i)*rtd enddo do i=1,nmlat gmlat(i) = ylatm(i)*rtd enddo do i=1,nlev pmlev(i+3) = plev(i) enddo do i=3,1,-1 pmlev(i) = pmlev(i+1)-dlev enddo ! ! Get login name: logname = ' ' call getenv('LOGNAME',logname) if (len_trim(logname)==0) then write(6,"(/,'>>> init: Cannot get LOGNAME environment ', | 'variable.',/)") stop 'LOGNAME' endif ! ! Get host name: call gethostsname(host) ! ! Operating system (based on pre-proc macro): call setosys(system) ! ! Get run date (current date and time): call datetime(rundate,runtime) ! ! Iter is the number of time steps from 0,0,0 to the current model ! time, using the current step length. Iter is incremented once per ! timestep in advnce. ! iter = mtime_to_nstep(start(:,1),step) ! ! iyear and iday are current calendar year and day. ! If calendar_advance > 0, the model is advanced in calendar time, ! starting at start_day. If calendar_advance==0, model is NOT ! advanced in calendar time (start_day is held constant). ! iyear and iday are incremented in advance if the model is ! advancing in calendar time. ! iyear = start_year ! from input iday = start_day ! from input sfeps = 1. ! ! If model is being advanced in calendar time, initialize orbital ! eccentricity. ! if (calendar_advance > 0) then theta0 = 2.*pi*float(iday)/365. sfeps = 1.000110+0.034221*cos(theta0)+0.001280*sin(theta0)+ | 0.000719*cos(2.*theta0)+0.000077*sin(2.*theta0) endif ! ! 2/00: these were in modsrc.snoe (tgcm13mt), but were unused. ! Low-energy protons: ! alfalp = 10. ! efluxlp = 1.e-20 ! ! Set GPI flag: igetgpi = 0 if (len_trim(gpi_ncfile) > 0) igetgpi = 1 if (igetgpi > 0) | write(6,"(' gpi_ncfile = ',a)") trim(gpi_ncfile) ! ! Set GSWM flag for diurnal tides: igetgswmdi = 0 if (len_trim(gswm_di_ncfile) > 0) igetgswmdi = 1 if (igetgswmdi > 0) | write(6,"(' gswm_di_ncfile = ',a)") trim(gswm_di_ncfile) ! Set GSWM flag for semidiurnal tides: igetgswmsdi = 0 if (len_trim(gswm_sdi_ncfile) > 0) igetgswmsdi = 1 if (igetgswmsdi > 0) | write(6,"(' gswm_sdi_ncfile = ',a)") trim(gswm_sdi_ncfile) ! Set GSWM flag for nonmigrating diurnal tides: igetgswmnmidi = 0 if (len_trim(gswm_nmidi_ncfile) > 0) igetgswmnmidi = 1 if (igetgswmnmidi > 0) | write(6,"(' gswm_nmidi_ncfile = ',a)") trim(gswm_nmidi_ncfile) ! Set GSWM flag for nonmigrating semidiurnal tides: igetgswmnmisdi = 0 if (len_trim(gswm_nmisdi_ncfile) > 0) igetgswmnmisdi = 1 if (igetgswmnmisdi > 0) | write(6,"(' gswm_nmisdi_ncfile = ',a)") | trim(gswm_nmisdi_ncfile) ! ! ixtimep is 4th dimension index to fg-array for previous time step ! ixtimec is 4th dimension index to fg-array for current time step ! (see fogcm.f) ! ! ixtimep = 1 ! ixtimec = 1 ! ! Initialize amie, and get amie file if necessary: ! call init_amie ! ! Initialize non-input history variables for beginning of run: call hist_init ! ! setfft calls set99 for fft init. This call returns trigs and ifax, ! in filter_module.F. ! call setfft(trigs,ifax,nlon) ! ! Init starting model time: start_mtime = start(:,1) ! ! Initialize field structures: iprintf = 0 call init_fields(lon0,lon1,lat0,lat1,mytid,iprintf) ! ! Initialize secondary history fields: if (isechist > 0) then call set_fsech endif ! ! Initialize sun's declination: sundec=atan(tan(23.5*pi/180.)*sin(2.*pi*float(iday-80)/365.)) sin_sundec = SIN(sundec) ! C(95) cos_sundec = COS(sundec) ! C(96) ! ! Define geographic grid: do i=1,nlon glon(i) = glon1+(i-1)*dlon enddo do i=1,nlat glat(i) = glat1+(i-1)*dlat enddo do i=1,nlev plev(i) = plev1+(i-1)*dlev enddo ! ! Report to stdout: write(6,"(/,'Model run initialization:')") write(6,"(' nstep = ',i6,4x, | '(Number of time steps this run)')") nstep write(6,"(' iter = ',i6,4x, | '(Initial iteration number)')") iter write(6,"(' iyear = ',i6,4x, | '(Beginning calendar year)')") iyear write(6,"(' iday = ',i6,4x, | '(Beginning calendar day)')") iday write(6,"(' igetgpi = ',i6,4x, | '(If > 0, geophysical indices database will be used.)')") | igetgpi write(6,"(' igetgswmdi = ',i6,4x, | '(If > 0, GSWM diurnal tidal database will be used.)')") | igetgswmdi write(6,"(' igetgswmsdi= ',i6,4x, | '(If > 0, GSWM semidiurnal tidal database will be used.)')") | igetgswmsdi write(6,"(' igetgswmnmidi= ',i6,4x, | '(If > 0, GSWM nonmigrating diurnal tidal database will', | ' be used.)')") igetgswmnmidi write(6,"(' igetgswmnmisdi= ',i6,4x, | '(If > 0, GSWM nonmigrating semidiurnal tidal database will', | ' be used.)')") igetgswmnmisdi ! ! ncep/nmc are in time-gcm only: ! write(6,"(' ncep = ',i6,4x, ! | '(If > 0, use NCEP Z and TN 10 mb lower boundaries.')") ! | ncep ! write(6,"(' nmc = ',i6,4x, ! | '(If > 0, use NMC Z and TN 10 mb lower boundaries.')") ! | nmc ! if (nsource > 0) then write(6,"(/,'This is an initial run:')") write(6,"(' start_year = ',i6,5x, | '(Starting year of initial run)')") start_year write(6,"(' start_day = ',i6,5x, | '(Starting day of initial run)')") start_day write(6,"(' start_mtime= ',i4,2i3,1x, | '(Starting mtime of initial run)')") start_mtime endif ! ! Report re primary histories to stdout: write(6,"(/,'Primary Histories:')") write(6,"(' nsource = ',i5,2x, | '(If > 0, a primary source history was provided)')") nsource write(6,"(' nseries_prim = ',i5,2x, | '(Number of primary time series)')") nseries_prim write(6,"(' nhist_total = ',i5,2x, | '(Number of primary histories to be written)')") nhist_total write(6,"(' nfiles_prim = ',i5,2x, | '(Number of primary output files to be written)')") nfiles_prim write(6,"(' mxhist_prim = ',i5,2x, | '(Maximum number of primary histories per file)')") mxhist_prim ! ! Report re secondary histories to stdout: if (isechist > 0) then write(6,"(/,'Secondary Histories:')") write(6,"(' nseries_sech = ',i5,2x, | '(Number of secondary time series)')") nseries_sech write(6,"(' nsech_total = ',i5,2x, | '(Number of secondary histories to be written)')") nsech_total write(6,"(' nfiles_sech = ',i5,2x, | '(Number of secondary output files to be written)')") | nfiles_sech write(6,"(' mxhist_sech = ',i5,2x, | '(Maximum number of secondary histories per file)')") | mxhist_sech write(6,"(' nfsech_geo = ',i5,2x, | '(Number of secondary history fields on geographic grid)')") | nfsech_geo write(6,"(' nfsech_mag = ',i5,2x, | '(Number of secondary history fields on magnetic grid)')") | nfsech_mag write(6,"(' nfsech_geo2d = ',i5,2x, | '(Number of secondary history fields on geographic ', | '2d grid)')")nfsech_geo2d write(6,"(' nfsech_mag2d = ',i5,2x, | '(Number of secondary history fields on magnetic 2d grid)')") | nfsech_mag2d write(6,"(' nfsech_magphr = ',i5,2x, | '(Number of secondary history fields on magnetoshere grid)')") | nfsech_magphr ! ! Report secondary history fields: write(6,"(/,'Secondary history fields on geographic grid ', | ' (number of fields =',i3,'):')") nfsech_geo do i=1,nfsech_geo if (len_trim(fsech(i)%long_name)>0) then write(6,"(' Field ',a,' (',a,')')") | fsech(i)%short_name(1:8),trim(fsech(i)%long_name) else write(6,"(' Field ',a)") fsech(i)%short_name(1:8) endif enddo write(6,"(/,'Secondary history fields on magnetic ', | 'grid (number of fields =',i3,'):')") nfsech_mag do i=1,nfsech_mag write(6,"(' Mag field ',a,' (diagnostic)')") | fsechmag(i)%short_name enddo write(6,"(/,'Secondary history fields on geographic 2d grid ', | ' (number of fields =',i3,'):')") nfsech_geo2d do i=1,nfsech_geo2d write(6,"(' Field ',a,' (diagnostic)')") | fsech2d(i)%short_name enddo write(6,"(/,'Secondary history fields on magnetic ', | '2d grid (number of fields =',i3,'):')") nfsech_mag2d do i=1,nfsech_mag2d write(6,"(' Mag 2d field ',a,' (diagnostic)')") | fsechmag2d(i)%short_name enddo write(6,"(/,'Secondary history fields on magnetospheric ', | '2d grid (number of fields =',i3,'):')") nfsech_magphr do i=1,nfsech_magphr write(6,"(' Magphr field ',a,' (diagnostic)')") | fsechmagphr2d(i)%short_name enddo endif end subroutine init !----------------------------------------------------------------------- subroutine init_geo2mag ! ! Determine points necessary for parallel geo to mag transformation. ! This is called from main tgcm after magfield (cannot be called from ! init). ! ! Module data from mpi.F: ! integer,allocatable,dimension(:,:,:) :: ! (mlon0->1,mlat0->1,4) ! | igeo, ! geo lon indices at 4 pts boxing each mag pt in my subdomain ! | jgeo, ! geo lat indices at 4 pts boxing each mag pt in my subdomain ! | idgeo ! tid of mpi task with each of the 4 pts (owner) ! integer,allocatable :: nrcvgeo(:) ! # geo pts to be recvd from each task ! integer,allocatable :: nsndgeo(:) ! # geo pts to be sent to each ntask ! integer,allocatable,dimension(:,:) :: ! (mxrcv[mxsnd],ntask) ! | ircvgeo, ! lon indices to be received from each task ! | jrcvgeo, ! lat indices to be received from each task ! | isndgeo, ! lon indices to be sent from each task ! | jsndgeo ! lat indices to be sent from each task ! use mpi_module,only: ntask,mlon0,mlon1,mlat0,mlat1,mymlons, | mymlats,igeo,jgeo,idgeo,ircvgeo,jrcvgeo,isndgeo,jsndgeo, | nsndgeo,nrcvgeo,tasks,ntaskj,mytid use magfield_module,only: ig, jg ! ! Local: integer :: im,jm,i,j,n,nmlats,nmlons,istat integer :: mxsndgeo,mxrcvgeo ! logical,dimension(nlonp4,nlat) :: counted logical,allocatable,dimension(:,:,:) :: | counted ! (nlonp4,nlat,ntask) ! ! Set magnetic grid in degrees: ! ! Allocate arrays (module data from mpi.F): allocate(igeo(mlon0:mlon1,mlat0:mlat1,4),stat=istat) if (istat /= 0) write(6,"('>>> init_geo2mag: error allocating', | ' igeo: stat=',i3)") istat allocate(jgeo(mlon0:mlon1,mlat0:mlat1,4),stat=istat) if (istat /= 0) write(6,"('>>> init_geo2mag: error allocating', | ' jgeo: stat=',i3)") istat allocate(idgeo(mlon0:mlon1,mlat0:mlat1,4),stat=istat) if (istat /= 0) write(6,"('>>> init_geo2mag: error allocating', | ' idgeo: stat=',i3)") istat allocate(nrcvgeo(0:ntask-1),stat=istat) if (istat /= 0) write(6,"('>>> init_geo2mag: error allocating', | ' nrcvgeo: stat=',i3)") istat allocate(nsndgeo(0:ntask-1),stat=istat) if (istat /= 0) write(6,"('>>> init_geo2mag: error allocating', | ' nsndgeo: stat=',i3)") istat ! ! Assign igeo and jgeo from ig,jg (which were read by magfield): ! Numbering convention for the 4 geo pts surrounding a single mag pt ! is 1 to 4, counterclockwise starting with 1 in lower left. ! do jm=mlat0,mlat1 do im=mlon0,mlon1 igeo(im,jm,1) = ig(im,jm) igeo(im,jm,2) = ig(im,jm)+1 igeo(im,jm,3) = ig(im,jm)+1 igeo(im,jm,4) = ig(im,jm) ! jgeo(im,jm,1) = jg(im,jm) jgeo(im,jm,2) = jg(im,jm) jgeo(im,jm,3) = jg(im,jm)+1 jgeo(im,jm,4) = jg(im,jm)+1 enddo ! im=mlon0,mlon1 enddo ! jm=mlat0,mlat1 ! nmlats = mlat1-mlat0+1 nmlons = mlon1-mlon0+1 ! ! Determine owners (tasks) of each mag grid pt in my subdomain: idgeo(:,:,:) = -1 do jm=mlat0,mlat1 do im=mlon0,mlon1 do i=1,4 do n=0,ntask-1 if (any(tasks(n)%mylons==igeo(im,jm,i)).and. | any(tasks(n)%mylats==jgeo(im,jm,i))) then idgeo(im,jm,i) = n ! this needed pt is owned by task n elseif (any(tasks(n)%mylons==igeo(im,jm,i)).and. | jgeo(im,jm,i) > nlat .and. | tasks(n)%mytidj==ntaskj-1) then ! northernmost tasks idgeo(im,jm,i) = n ! owner of nlat+1,nlat+2 endif enddo ! n=0,ntask-1 if (idgeo(im,jm,i)==-1) then write(6,"('>>> init_geo2mag note: no owner: ', | ' jm=',i3,' im=',i3,' i=',i3,' igeo=',i3,' jgeo=',i3)") | jm,im,i,igeo(im,jm,i),jgeo(im,jm,i) stop 'init_geo2mag' endif enddo ! i=1,4 enddo ! im=mlon0,mlon1 enddo ! jm=mlat0,mlat1 if (any(idgeo==-1)) then write(6,"('>>> WARNING init_geo2mag: some idgeo == -1')") stop 'idgeo' endif ! ! Determine how many points I must receive from each task, for ! allocation purposes below. ! Use logical counted to remove redundancy -- only need to send/recv ! each point once. ! allocate(counted(nlonp4,nlat,0:ntask-1),stat=istat) if (istat /= 0) write(6,"('>>> init_geo2mag: error allocating', | ' counted: stat=',i3)") istat counted(:,:,:) = .false. ! (nlonp4,nlat,ntask) nrcvgeo(:) = 0 nsndgeo(:) = 0 do jm=mlat0,mlat1 do im=mlon0,mlon1 do i=1,4 n = idgeo(im,jm,i) if (.not.counted(ig(im,jm),jg(im,jm),n)) then nrcvgeo(n) = nrcvgeo(n)+1 if (n==mytid) nsndgeo(mytid) = nsndgeo(mytid)+1 counted(ig(im,jm),jg(im,jm),n) = .true. endif enddo ! i=1,4 enddo ! im=mlon0,mlon1 enddo ! jm=mlat0,mlat1 do n=0,ntask-1 write(6,"('init_geo2mag: n=',i3,' nrcvgeo=',i4,' nsndgeo=',i4)") | n,nrcvgeo(n),nsndgeo(n) enddo ! ! Allocate ircvgeo,jrcvgeo: mxrcvgeo = 0 do n=0,ntask-1 if (nrcvgeo(n) > mxrcvgeo) mxrcvgeo = nrcvgeo(n) if (nsndgeo(n) > mxsndgeo) mxsndgeo = nsndgeo(n) enddo allocate(ircvgeo(mxrcvgeo,0:ntask-1),stat=istat) if (istat /= 0) write(6,"('>>> init_geo2mag: error allocating', | ' ircvgeo: stat=',i3)") istat allocate(jrcvgeo(mxrcvgeo,0:ntask-1),stat=istat) if (istat /= 0) write(6,"('>>> init_geo2mag: error allocating', | ' jrcvgeo: stat=',i3)") istat allocate(isndgeo(mxsndgeo,0:ntask-1),stat=istat) if (istat /= 0) write(6,"('>>> init_geo2mag: error allocating', | ' isndgeo: stat=',i3)") istat allocate(jsndgeo(mxsndgeo,0:ntask-1),stat=istat) if (istat /= 0) write(6,"('>>> init_geo2mag: error allocating', | ' jsndgeo: stat=',i3)") istat ! ! Determine lon,lat indices to be passed by each task. ! (recalc nrcvgeo, nsndgeo) counted(:,:,:) = .false. ! (nlonp4,nlat,ntask) nrcvgeo(:) = 0 nsndgeo(:) = 0 write(6,"('needed fm',2x,'nrcvgeo',2x,' mag lon',2x, | ' mag lat',2x,' geo lon',2x,' geo lat')") ! 11/21/02: Hard to prove that this is really working.. do jm=mlat0,mlat1 do im=mlon0,mlon1 do i=1,4 n = idgeo(im,jm,i) if (.not.counted(ig(im,jm),jg(im,jm),n)) then nrcvgeo(n) = nrcvgeo(n)+1 ircvgeo(nrcvgeo(n),n) = igeo(im,jm,i) jrcvgeo(nrcvgeo(n),n) = jgeo(im,jm,i) if (n==mytid) then nsndgeo(mytid) = nsndgeo(mytid)+1 isndgeo(nsndgeo(n),n) = igeo(im,jm,i) jsndgeo(nsndgeo(n),n) = jgeo(im,jm,i) endif counted(ig(im,jm),jg(im,jm),n) = .true. write(6,"(i9,2x,i7,2x,f9.2,2x,f9.2,2x,f9.2,2x, | f9.2)") n,nrcvgeo(n),glon(ircvgeo(nrcvgeo(n),n)), | glat(ircvgeo(nrcvgeo(n),n)),gmlon(im),gmlat(jm) endif enddo ! i=1,4 enddo ! im=mlon0,mlon1 enddo ! jm=mlat0,mlat1 deallocate(counted) ! end subroutine init_geo2mag !----------------------------------------------------------------------- end module init_module