!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> Module used to exchange data via InterComm with coupled models !! - MIX (Coupler) !! - LFM (MHD) !! !! FIXME: Unit conversion in the Import & Export functions should be !! moved to some generic prepareExport & processImport function. The !! intercomm_module should only contain code for communication. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! module intercomm_module use dipole_params, only : dm !> Dipole Moment use lfm_transfergrid, only : setupLFMGrid, tearDownLFMGrid, & ni_mid, nj_mid, nk_mid, & xmid_lfm, ymid_lfm, zmid_lfm, & bleed_RHO_lfm, & !> RCM's Bleed Density on LFM grid bleed_PRESSURE_lfm, & !> RCM's Bleed Pressure on LFM grid bleed_MASK_lfm !> RCM's Bleed Mask on LFM grid: 0=not part of RCM; 1=part of RCM use intermediate_grid, only : idim_ig, jdim_ig, kdim_ig, & x_ig, y_ig, z_ig, & bx_ig, by_ig, bz_ig, & !> Accumulate magnetic field press_ig, & !> Accumulate pressure rho_ig, & !> Accumulate density vel_ig !> Accumulate speed use ionosphere_intermediate, only : nLat_ion, nLon_ion, & gcolat, glong, & pot, & !> Potential received from MIX coupler sigmap, & !> Pedersen Conductance received from MIX coupler sigmah, & !> Hall Conductance received from MIX coupler eng_avg, & !> Average Energy sent to MIX coupler flux, & !> Average energy flux sent to MIX coupler fac !> Total field-aligned current density sent to MIX coupler implicit none !> Public functions/subroutines public Initialize, ReceiveScalars, & ExportIonosphere, ExportMagnetosphere, & ImportIonosphere, ImportMagnetosphere, & Finalize !> Date & Time variables integer, public :: iaUtRcmStart(6) !> date/time to start coupling with MHD integer, public :: iRcmDelta !> Frequency (in seconds) to exchange data with MHD & CPL !> Indices for status scalars to be sent from the MHD !! Note: These must match up with the scalars on the sending side. integer, parameter, public :: KILL_SIGNAL = 1 !> 111=kil; 222=last exchange; else=keep running integer, parameter, public :: YEAR = 2 integer, parameter, public :: MONTH = 3 integer, parameter, public :: DAY = 4 integer, parameter, public :: HOUR = 5 integer, parameter, public :: MINUTE = 6 integer, parameter, public :: SECOND = 7 integer, parameter, public :: DELTA_T = 8 !> # of seconds between exchange with coupled models integer, parameter, public :: LABEL = 9 !> MHD Simulation Time step integer, parameter, public :: NUMBER_OF_SCALARS = 9 ! Scalars at index KILL_SIGNAL accept the following: integer, parameter, public :: KILL_SIGNAL_CONTINUE = 0 integer, parameter, public :: KILL_SIGNAL_SHUTDOWN = 111 integer, parameter, public :: KILL_SIGNAL_LAST_EXCHANGE = 222 private !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> Constants character(LEN=3), parameter :: prog_name = "rcm" !> program name which must match XJD file value real :: fRi !> FIXME: What is this, again? !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> Intermediate Grid accumulate (from LFM to RCM) variables !! unwrapped as 1d arrays. Internal RCM arrays defined in the !! "intermediate_grid" module. ! FIXME: Should refactor code to remove these. Grep this src for "RESHAPE" for details. real, allocatable :: bx_ig_1d(:), by_ig_1d(:), bz_ig_1d(:), & press_ig_1d(:), rho_ig_1d(:), vel_ig_1d(:) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> InterComm variables integer :: xjd !> xjd id value integer :: status !> return status from Intercomm calls ! InterComm descriptors integer :: descriptor_mhd !> lfm {x,y,z} grid, density, pressure, mask integer :: descriptor_igxyz !> accumulate density, pressure, b{x,y,z}, v integer :: descriptor_lat !> RCM ionosphere grid latitudes integer :: descriptor_lon !> RCM ionosphere grid longitudes integer :: descriptor_ion !> avg_energy, avg_flux, fac, potential, conductances ! InterComm Regions integer :: regions_mhd !> lfm {x,y,z} grid, density, pressure, mask integer :: regions_igxyz !> accumulate density, pressure, b{x,y,z}, v integer :: regions_lat !> RCM ionosphere grid latitudes integer :: regions_lon !> RCM ionosphere grid longitudes integer :: regions_ion !> avg_energy, avg_flux, fac, potential, conductances contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> Initialize !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine Initialize use constants, only : radius_earth_m, & !> used to scale dipole moment radius_iono_m use rice_housekeeping_module, only : xjd_filename !> XJD file that describes all connections implicit none !> Constants integer, parameter :: IC_COLUMN_MAJOR = 1 !> Shouldn't this be declared somewhere in InterComm? !> Helper variables... integer :: i,j,k !> Loop indices integer :: iaScalars(NUMBER_OF_SCALARS) !> Status scalars sent from MHD integer :: nIgPoints !> Number of intermediate grid points (idim_ig * jdim_ig * kdim_ig) integer :: iaIgSize(3) = (/-1, -1, -1 /) !> Intermediate Grid dimensions {idim_ig,jdim_ig,kdim_ig} integer :: iaRcmIonDimensions (2) !> Dimensions of RCM ionosphere lat/lon grid integer :: iaLFMSize(3)=(/ -1, -1, -1/) !> LFM grid dimensions (ni_mid, nj_mid, nk_mid) !> InterComm Variables... integer :: rank=0 !> since RCM is serial, this is 0 for now integer :: iTasks(1) !> array holding assignments of domain decomposition integer :: iaBlocks1d(1,2,1), iaLower1d(1), iaUpper1d(1), iaStride1d(1) !> 1d regions integer :: iaBlocks2d(1,2,2), iaLower2d(2), iaUpper2d(2), iaStride2d(2) !> 2d regions integer :: iaBlocks3d(1,2,3), iaLower3d(3), iaUpper3d(3), iaStride3d(3) !> 3d regions !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! CALL IC_Initialize (prog_name, rank, xjd_filename, xjd, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_Initialize") ! This is a serial code, so tasks=0 (note: according to InterComm manual, we start from 0, even in Fortran.) iTasks (1) = 0 ! Send Intermediate Grid dimensinos. . . nIgPoints=idim_ig*jdim_ig*kdim_ig iaIgSize = (/ idim_ig, jdim_ig, kdim_ig /) CALL IC_bcast_local (xjd, "im-mhd_ig_size", iaIgSize, 3, status) ! send IG sizes to LFM CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_bcast_local") ! Send RCM Ionosphere Grid dimensions iaRcmIonDimensions = (/ nLat_ion, nLon_ion /) CALL IC_bcast_local (xjd, "im-cpl_im_size", iaRcmIonDimensions, 2, status) ! send RCM grid sizes to LFM CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_bcast_local") CALL IC_Recv_local (xjd, "mhd-im_lfm_size", iaLFMSize, 3, status) ! get LFM grid sizes from LFM CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_recv_local") print*,"RCM received LFM grid sizes ni_mid,nj_mid, nk_mid=: ", iaLFMSize ! Arrays must be created at this point, before defining regions ! with Intercomm call setupLFMGrid(iaLFMSize(1), iaLFMSize(2), iaLFMSize(3)) ALLOCATE (rho_ig_1d(nIgPoints)) ALLOCATE (press_ig_1d(nIgPoints)) ALLOCATE (bx_ig_1d (nIgPoints)) ALLOCATE (by_ig_1d (nIgPoints)) ALLOCATE (bz_ig_1d (nIgPoints)) ALLOCATE (vel_ig_1d (nIgPoints)) write(6,*) "Allocated memory. Continue exchanging variables" ! Send/receive scalar variables: CALL IC_Recv_local (xjd, "mhd-im_dipole_moment", dm, 1, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_recv_local") print*, "RCM received dipole moment=", dm ! convert Dipole Moment to nt/re^3 !FIXME: Check units dm = dm / (radius_earth_m**3)*1.0e+9 ! Receive radius of ionosphere from MIX & make sure RCM has been ! using the same value. Note: RCM requires rIon before InterComm ! is initialized, so this is just some error checking... fRi = 0. CALL IC_Recv_local (xjd, "cpl-im_r_i", fRi, 1, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_bcast_local") IF ( abs(fRi-radius_iono_m) .gt. 1e-6 ) then STOP " intercomm.f90: Radius of ionosphere differs between RCM and the Coupler (MIX). These values must be identical!" END IF write(*,*) "Broadcasting ", idim_ig, jdim_ig, kdim_ig ! Convert Intermediate Grid from Earth Radii to MKS ! for the purpose of sending it to LFM only: x_ig = x_ig * radius_earth_m y_ig = y_ig * radius_earth_m z_ig = z_ig * radius_earth_m ! Send intermediate grid to LFM: CALL IC_bcast_local (xjd, "im-mhd_igx", x_ig, idim_ig, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_Export(IGX)") CALL IC_bcast_local (xjd, "im-mhd_igy", y_ig, jdim_ig, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_Export(IGY)") CALL IC_bcast_local (xjd, "im-mhd_igz", z_ig, kdim_ig, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_Export(IGZ)") ! Convert Intermediate Grid from MKS back into Earth Radii (native units for RCM) ! after it has been sent to LFM: x_ig = x_ig / radius_earth_m y_ig = y_ig / radius_earth_m z_ig = z_ig / radius_earth_m write(6,*) "Done exchanging initial variables." ! Register 3-D LFM arrays with Intercomm: iaBlocks3d(1,1:2,1) = (/1, ni_mid /) iaBlocks3d(1,1:2,2) = (/1, nj_mid /) iaBlocks3d(1,1:2,3) = (/1, nk_mid /) iaLower3d = (/1,1,1 /) iaUpper3d = (/ ni_mid, nj_mid, nk_mid /) iaStride3d = (/1,1,1/) CALL IC_Create_bdecomp_desc (3, iaBlocks3d, iTasks, 1, descriptor_mhd, IC_COLUMN_MAJOR) CALL IC_Create_block_region (3, iaLower3d, iaUpper3d, iaStride3d, regions_mhd) CALL IC_Register_region (xjd, regions_mhd, 1, "mhd-im_lfmx", descriptor_mhd, xmid_lfm, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(LFMX)") CALL IC_Register_region (xjd, regions_mhd, 1, "mhd-im_lfmy", descriptor_mhd, ymid_lfm, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(LFMY)") CALL IC_Register_region (xjd, regions_mhd, 1, "mhd-im_lfmz", descriptor_mhd, zmid_lfm, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(LFMZ)") CALL IC_Register_region (xjd, regions_mhd, 1, "im-mhd_density", descriptor_mhd, bleed_RHO_lfm, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(BLEEDRHO)") CALL IC_Register_region (xjd, regions_mhd, 1, "im-mhd_pressure", descriptor_mhd, bleed_PRESSURE_lfm, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(BLEEDP)") CALL IC_Register_region (xjd, regions_mhd, 1, "im-mhd_rcm_mask", descriptor_mhd, bleed_MASK_lfm, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(BLEEDP)") ! Register Intermediate Grid arrays (packed into 1-D for LFM): iaBlocks1d(1,1:2,1) = (/ 1, nIgPoints /) iaLower1d (1) = 1 iaUpper1d (1) = nIgPoints iaStride1d(1) = 1 CALL IC_Create_bdecomp_desc (1, iaBlocks1d, iTasks, 1, descriptor_igxyz, IC_COLUMN_MAJOR) CALL IC_Create_block_region (1, iaLower1d, iaUpper1D, iaStride1D, regions_igxyz) CALL IC_Register_region (xjd, regions_igxyz, 1, "mhd-im_density", descriptor_igxyz, rho_ig_1d, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(ACCUMRHO)") CALL IC_Register_region (xjd, regions_igxyz, 1, "mhd-im_pressure", descriptor_igxyz, press_ig_1d, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(ACCUMP)") CALL IC_Register_region (xjd, regions_igxyz, 1, "mhd-im_bx", descriptor_igxyz, bx_ig_1d, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(ACCUMBX)") CALL IC_Register_region (xjd, regions_igxyz, 1, "mhd-im_by", descriptor_igxyz, by_ig_1d, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(ACCUMBY)") CALL IC_Register_region (xjd, regions_igxyz, 1, "mhd-im_bz", descriptor_igxyz, bz_ig_1d, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(ACCUMBZ)") CALL IC_Register_region (xjd, regions_igxyz, 1, "mhd-im_v", descriptor_igxyz, vel_ig_1d, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(ACCUMV)") ! Register RCM ionospheric grid arrays: iaBlocks1d(1,1:2,1) = (/ 1, nLat_ion /) iaLower1d (1) = 1 iaUpper1d (1) = nLat_ion iaStride1d(1) = 1 CALL IC_Create_bdecomp_desc (1, iaBlocks1d, iTasks, 1, descriptor_lat, IC_COLUMN_MAJOR) CALL IC_Create_block_region (1, iaLower1d, iaUpper1D, iaStride1D, regions_lat) CALL IC_Register_region (xjd, regions_lat, 1, "im-cpl_lat", descriptor_lat, gcolat, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(RCMIONLAT)") iaBlocks1d(1,1:2,1) = (/ 1, nLon_ion /) iaLower1d (1) = 1 iaUpper1d (1) = nLon_ion iaStride1d(1) = 1 CALL IC_Create_bdecomp_desc (1, iaBlocks1d, iTasks, 1, descriptor_lon, IC_COLUMN_MAJOR) CALL IC_Create_block_region (1, iaLower1d, iaUpper1D, iaStride1D, regions_lon) CALL IC_Register_region (xjd, regions_lon, 1, "im-cpl_lon", descriptor_lon, glong, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_register_region(RCMIonLong)") iaBlocks2d (1,1:2,1) = (/ 1, nLat_ion /) iaBlocks2d (1,1:2,2) = (/ 1, nLon_ion /) iaLower2d (:) = (/ 1, 1 /) iaUpper2d (:) = (/ nLat_ion, nLon_ion /) iaStride2d = (/1,1/) CALL IC_Create_bdecomp_desc (2, iaBlocks2d, iTasks, 1, descriptor_ion, IC_COLUMN_MAJOR) CALL IC_Create_block_region (2,iaLower2d,iaUpper2d, iaStride2d, regions_ion) CALL IC_Register_region (xjd, regions_ion, 1, "im-cpl_avg_energy", descriptor_ion, eng_avg, status) CALL CheckIcStatus (xjd, status, " RCM_MIX/IC_register_region(avgEnergy)") CALL IC_Register_region (xjd, regions_ion, 1, "im-cpl_avg_flux", descriptor_ion, flux, status) CALL CheckIcStatus (xjd, status, " RCM_MIX/IC_register_region(avgFlux)") CALL IC_Register_region (xjd, regions_ion, 1, "im-cpl_fac", descriptor_ion, fac, status) CALL CheckIcStatus (xjd, status, " RCM_MIX/IC_register_region(avgFAC)") CALL IC_Register_region (xjd, regions_ion, 1, "cpl-im_potential", descriptor_ion, pot, status) CALL CheckIcStatus (xjd, status, " MIX_RCM/IC_register_region(Potential)") CALL IC_Register_region (xjd, regions_ion, 1, "cpl-im_sigma_p", descriptor_ion, sigmap, status) CALL CheckIcStatus (xjd, status, " MIX_RCM/IC_register_region(sigma_p)") CALL IC_Register_region (xjd, regions_ion, 1, "cpl-im_sigma_h", descriptor_ion, sigmah, status) CALL CheckIcStatus (xjd, status, " MIX_RCM/IC_register_region(sigma_h)") write(6,*) "Regions defined. Comitting." ! Finish registering arrays by committing them: CALL IC_Commit_region (xjd, status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_Commit_region") !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! write(6,*) "Comitted Arrays. Handshaking." !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Send RCM ionospheric grid to LFM: CALL IC_Export (xjd, "im-cpl_lat", status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_Export(RCMIONLAT)") CALL IC_Export (xjd, "im-cpl_lon", status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_Export(RCMIONLONG)") ! Receive from LFM its grid: CALL IC_Import (xjd, "mhd-im_lfmx", status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_Import(LFMX)") CALL IC_Import (xjd, "mhd-im_lfmy", status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_Import(LFMY)") CALL IC_Import (xjd, "mhd-im_lfmz", status) CALL CheckIcStatus (xjd, status, " RCM_LFM/IC_Import(LFMZ)") ! Convert LFM grid midpoints to units of Re as needed by RCM: xmid_lfm = xmid_lfm / radius_earth_m ymid_lfm = ymid_lfm / radius_earth_m zmid_lfm = zmid_lfm / radius_earth_m ! ___ Debugging output _________________________________________ ! open(unit=65,file='lfm_midpoint_grid.dat',status='unknown') ! ! Dimensions: ! write(65,*) ni_mid, nj_mid, nk_mid ! write(65,*) "i j k lfmx lfmy lfmz" ! do i=1, ni_mid ! do j=1, nj_mid ! do k=1, nk_mid ! write(65,*) i-1, j-1, k-1, & ! XMID_lfm(i,j,k), & ! YMID_lfm(i,j,k), & ! ZMID_lfm(i,j,k) ! end do ! end do ! end do ! close(65) ! ___ Debugging output _________________________________________ ! open(unit=65,file='GRIDOUT_rcm.dat',status='unknown') ! ! Dimensions: ! write(65,*) ni_mid, nj_mid, -999, -999, nk_mid ! ! PHI: ! do k=1, nk_mid ! write(65,*) atan2(zmid_lfm(10,10,k), ymid_lfm(10,10,k))*180.0/3.14159265 ! end do ! ! X2: ! do j=1, nj_mid ! do i=1, ni_mid ! write(65,*) xmid_lfm(i,j,1) ! end do ! end do ! ! Y2: ! do j=1, nj_mid ! do i=1, ni_mid ! write(65,*) sqrt( ymid_lfm(i,j,1)**2 + zmid_lfm(i,j,1)**2 ) ! end do ! end do ! close(65) ! Receive from LFM RCM start time and Delta_T: call ReceiveScalars("mhd-im_startup_scalars", iaScalars) iaUtRcmStart(1) = iaScalars(YEAR) iaUtRcmStart(2) = iaScalars(MONTH) iaUtRcmStart(3) = iaScalars(DAY) iaUtRcmStart(4) = iaScalars(HOUR) iaUtRcmStart(5) = iaScalars(MINUTE) iaUtRcmStart(6) = iaScalars(SECOND) iRcmDelta = iaScalars(DELTA_T) end subroutine Initialize !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> Receive Scalars !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine ReceiveScalars(name, scalars) implicit none character (len = *), intent(in) :: name !> "mhd-im_scalars" or "mhd-im_startup_scalars" integer, intent(out) :: scalars(NUMBER_OF_SCALARS) CALL IC_Recv_local (xjd, name, scalars, NUMBER_OF_SCALARS, status) CALL CheckIcStatus (xjd, status, " IC_Recv_local(scalars)") end subroutine ReceiveScalars !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> Export data to Ionosphere Coupling (MIX) - Energy, Flux & FAC !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine ExportIonosphere ! ___ Unit conversion _________________________________________ ! FIXME: Are units in MKS? ! ___ Export data ______________________________________________ call print_max_2d(nLat_ion,nLon_ion,' Avg Energy',eng_avg,gcolat,glong) call print_min_2d(nLat_ion,nLon_ion,' Avg Energy',eng_avg,gcolat,glong) call print_max_2d(nLat_ion,nLon_ion,' Energy Flux',flux,gcolat,glong) call print_min_2d(nLat_ion,nLon_ion,' Energy Flux',flux,gcolat,glong) call print_max_2d(nLat_ion,nLon_ion,' FAC',fac,gcolat,glong) call print_min_2d(nLat_ion,nLon_ion,' FAC',fac,gcolat,glong) CALL IC_Export (xjd, "im-cpl_avg_energy", status) CALL IC_Export (xjd, "im-cpl_avg_flux", status) CALL IC_Export (xjd, "im-cpl_fac", status) end subroutine ExportIonosphere !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> Export data to Magnetosphere/MHD (LFM) - Bleed variables !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine ExportMagnetosphere integer :: i,j,k ! ___ Unit conversion _________________________________________ ! FIXME: Are units in MKS? ! ___ Export data ______________________________________________ call print_max(ni_mid,nj_mid,nk_mid,' im_density (kg/m^3)',bleed_RHO_lfm,XMID_lfm,YMID_lfm,ZMID_lfm) call print_min(ni_mid,nj_mid,nk_mid,' im_density (kg/m^3)',bleed_RHO_lfm,XMID_lfm,YMID_lfm,ZMID_lfm) call print_max(ni_mid,nj_mid,nk_mid,' im_pressure (Pa)',bleed_PRESSURE_lfm,XMID_lfm,YMID_lfm,ZMID_lfm) call print_min(ni_mid,nj_mid,nk_mid,' im_pressure (Pa)',bleed_PRESSURE_lfm,XMID_lfm,YMID_lfm,ZMID_lfm) call print_max_integer(ni_mid,nj_mid,nk_mid,' rcm mask (boolean)',bleed_MASK_lfm,XMID_lfm,YMID_lfm,ZMID_lfm) call print_min_integer(ni_mid,nj_mid,nk_mid,' rcm mask (boolean)',bleed_MASK_lfm,XMID_lfm,YMID_lfm,ZMID_lfm) CALL IC_Export (xjd, "im-mhd_density", status) CALL IC_Export (xjd, "im-mhd_pressure", status) CALL IC_Export (xjd, "im-mhd_rcm_mask", status) ! ___ Debugging output _________________________________________ ! open(unit=65,file='bleed_rcm.dat',status='unknown') ! write(65,*) ni_mid, nj_mid, nk_mid ! write(65,*) "i j k lfmx lfmy lfmz im_density im_pressure im_mask" ! do i=1, ni_mid ! do j=1, nj_mid ! do k=1, nk_mid ! write(65,*) i-1, j-1, k-1, & ! XMID_lfm(i,j,k), & ! YMID_lfm(i,j,k), & ! ZMID_lfm(i,j,k), & ! bleed_RHO_lfm(i,j,k), & ! bleed_PRESSURE_lfm(i,j,k), & ! bleed_MASK_lfm(i,j,k) ! end do ! end do ! end do ! close(65) end subroutine ExportMagnetosphere !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> Import data from Ionosphere Coupler (MIX) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine ImportIonosphere integer :: i,j real :: x,y ! ___ Import data ______________________________________________ CALL IC_Import (xjd, "cpl-im_potential", status) CALL IC_Import (xjd, "cpl-im_sigma_p", status) CALL IC_Import (xjd, "cpl-im_sigma_h", status) call print_max_2d(nLat_ion,nLon_ion,' Potential',pot,gcolat,glong) call print_min_2d(nLat_ion,nLon_ion,' Potential',pot,gcolat,glong) call print_max_2d(nLat_ion,nLon_ion,' Pedersen',sigmap,gcolat,glong) call print_min_2d(nLat_ion,nLon_ion,' Pedersen',sigmap,gcolat,glong) call print_max_2d(nLat_ion,nLon_ion,' Hall',sigmah,gcolat,glong) call print_min_2d(nLat_ion,nLon_ion,' Hall',sigmah,gcolat,glong) ! ___ Unit conversion _________________________________________ ! FIXME: Are units in MKS? ! ___ Debugging output _________________________________________ ! open(unit=65,file='rcm_rcm_grid.dat',status='unknown') ! write(65,*) nlon_ion, nlat_ion ! write(65,*) "i j x y potential sigma_p sigma_h" ! do i=1, nlon_ion ! do j=1, nlat_ion ! write(65,*) i-1, j-1, & ! sin(gcolat(j)) * cos(glong(i)), & ! sin(gcolat(j)) * sin(glong(i)), & ! pot(j,i), & ! sigmap(j,i), & ! sigmah(j,i) ! end do ! end do ! close(65) end subroutine ImportIonosphere !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> Import data from Magnetosphere/MHD (LFM) - Accumulate variables !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine ImportMagnetosphere ! ___ Import data ______________________________________________ CALL IC_Import (xjd, "mhd-im_density", status) CALL IC_Import (xjd, "mhd-im_pressure", status) CALL IC_Import (xjd, "mhd-im_bx", status) CALL IC_Import (xjd, "mhd-im_by", status) CALL IC_Import (xjd, "mhd-im_bz", status) CALL IC_Import (xjd, "mhd-im_v", status) ! ___ Unit conversion _________________________________________ ! density (rho_ig) is already in kg/m^3, no unit conversion. ! pressure (press_ig) is already in Pa, no unit conversion. ! Convert magnetic field from MKS units (Tesla) to nT bx_ig_1d = bx_ig_1d*1.0e9 by_ig_1d = by_ig_1d*1.0e9 bz_ig_1d = bz_ig_1d*1.0e9 ! Convert velocity to km/s vel_ig_1d = vel_ig_1d /1.0e3 ! in km/s ! ___ Unpack Arrays ___________________________________________ ! Unpack 1d intermediate grid arrays into 3 dimensions, since ! that's how RCM deals with IG vars internally. ! FIXME: This seems like a waste of memory... Refactor the ! code to use 1d arrays internally. press_ig = RESHAPE( press_ig_1d, (/idim_ig, jdim_ig, kdim_ig /) ) rho_ig = RESHAPE( rho_ig_1d, (/idim_ig, jdim_ig, kdim_ig /) ) bx_ig = RESHAPE( bx_ig_1d, (/idim_ig, jdim_ig, kdim_ig /) ) by_ig = RESHAPE( by_ig_1d, (/idim_ig, jdim_ig, kdim_ig /) ) bz_ig = RESHAPE( bz_ig_1d, (/idim_ig, jdim_ig, kdim_ig /) ) vel_ig = RESHAPE( vel_ig_1d, (/idim_ig, jdim_ig, kdim_ig /) ) end subroutine ImportMagnetosphere !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> Finalize: Free dynamically-allocated memory. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine Finalize ! Deallocate arrays call tearDownLFMGrid if (ALLOCATED(rho_ig_1d)) DEALLOCATE (rho_ig_1d) if (ALLOCATED(press_ig_1d)) DEALLOCATE (press_ig_1d) if (ALLOCATED(bx_ig_1d)) DEALLOCATE (bx_ig_1d) if (ALLOCATED(by_ig_1d)) DEALLOCATE (by_ig_1d) if (ALLOCATED(bz_ig_1d)) DEALLOCATE (bz_ig_1d) if (ALLOCATED(vel_ig_1d)) DEALLOCATE (vel_ig_1d) call IC_Free_desc(descriptor_mhd) call IC_Free_desc(descriptor_igxyz) call IC_Free_desc(descriptor_lat) call IC_Free_desc(descriptor_lon) call IC_Free_desc(descriptor_ion) call IC_Free_region(regions_mhd) call IC_Free_region(regions_igxyz) call IC_Free_region(regions_lat) call IC_Free_region(regions_lon) call IC_Free_region(regions_ion) call IC_Finalize (xjd, status) call CheckIcStatus(xjd, status, "Trouble with IC_Finalize") end subroutine Finalize !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> CheckIcStatus: Look at InterComm's error code & die if there !! are errors. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! SUBROUTINE CheckIcStatus (xjd,status, prefix_char) IMPLICIT NONE INTEGER, INTENT (IN) :: xjd,status CHARACTER (LEN=*), INTENT (IN) :: prefix_char INTEGER :: ist IF (status /= 0) then WRITE (*,*) prefix_char WRITE (*,'(A,I5)') " Status code: ",status WRITE (*,*) " ABORTING RCM_LFM" CALL IC_Finalize (xjd, ist) STOP END IF RETURN END SUBROUTINE CheckIcStatus end module intercomm_module