! $Id$ ! ! Earth System Modeling Framework ! Copyright (c) 2002-2023, University Corporation for Atmospheric Research, ! Massachusetts Institute of Technology, Geophysical Fluid Dynamics ! Laboratory, University of Michigan, National Centers for Environmental ! Prediction, Los Alamos National Laboratory, Argonne National Laboratory, ! NASA Goddard Space Flight Center. ! Licensed under the University of Illinois-NCSA License. ! !============================================================================== ! program ESMF_MeshOpUTest !------------------------------------------------------------------------------ #include "ESMF.h" !============================================================================== !BOPI ! !PROGRAM: ESMF_MeshOpUTest - Unit tests for Field Create and Get methods ! ! !DESCRIPTION: ! ! The code in this file drives F90 Field Create and Get unit tests. ! The companion folder Fieldsrc contains the definitions for the ! Field methods. !EOPI !----------------------------------------------------------------------------- ! !USES: use ESMF_TestMod ! test methods use ESMF implicit none !------------------------------------------------------------------------------ ! The following line turns the CVS identifier string into a printable variable. character(*), parameter :: version = & '$Id$' ! cumulative result: count failures; no failures equals "all pass" integer :: result = 0 ! individual test result code integer :: rc = 1 ! individual test failure message character(ESMF_MAXSTR) :: failMsg character(512) :: name call ESMF_TestStart(ESMF_SRCLINE, rc=rc) if (rc /= ESMF_SUCCESS) call ESMF_Finalize(endflag=ESMF_END_ABORT) !------------------------------------------------------------------------ !NEX_UTest ! exercise difference operation print *, 'Starting test3' call test3(rc) write(failMsg, *) "" write(name, *) "exercise difference operation" call ESMF_Test((rc.eq.ESMF_SUCCESS), name, failMsg, result, ESMF_SRCLINE) call ESMF_TestEnd(ESMF_SRCLINE) contains #define ESMF_METHOD "ESMF_TESTS" !------------------------------------------------------------------------ subroutine test3(rc) integer, intent(out) :: rc integer :: localrc type(ESMF_Mesh) :: meshA, meshB, mesh rc = ESMF_SUCCESS localrc = ESMF_SUCCESS meshA = make_mesh_sph(40,20,1.,1.,30.,-5.,rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return meshB = make_mesh_sph(20,10,0.5,1.,50.,0.,rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return mesh = ESMF_MeshCreate(meshA, meshB, ESMF_MESHOP_DIFFERENCE, & areaThreshold=1.e-7_ESMF_KIND_R8, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return !call ESMF_MeshWrite(mesh, "meshdiff", rc=localrc) !if (ESMF_LogFoundError(localrc, & ! ESMF_ERR_PASSTHRU, & ! ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_MeshDestroy(mesh, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return end subroutine test3 !------------------------------------------------------------------------ ! Utility methods !------------------------------------------------------------------------ function make_grid(atm_nx, atm_ny, atm_dx, atm_dy, atm_sx, atm_sy, field, rc) ! return value type(ESMF_Grid) :: make_grid ! arguments integer, intent(in) :: atm_nx, atm_ny real(ESMF_KIND_R4), intent(in) :: atm_dx, atm_dy real(ESMF_KIND_R4), intent(in) :: atm_sx, atm_sy type(ESMF_Field), intent(out), optional :: field integer, intent(out), optional :: rc ! local variables integer :: localrc, i, j real(ESMF_KIND_R8), pointer :: coordX(:), coordY(:) real(ESMF_KIND_R8) :: startx, starty make_grid = ESMF_GridCreateNoPeriDim(maxIndex=(/atm_nx, atm_ny/), & coordSys=ESMF_COORDSYS_CART, & indexflag=ESMF_INDEX_GLOBAL, & coordDep1=(/1/), & coordDep2=(/2/), & rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_GridAddCoord(make_grid, staggerloc=ESMF_STAGGERLOC_CENTER, & rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_GridAddCoord(make_grid, staggerloc=ESMF_STAGGERLOC_CORNER, & rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! global indexing ! atm grid is not decomposed in the y direction !startx = lpet*atm_nx/npet*atm_dx startx = atm_sx starty = atm_sy ! compute coord ! X center !call ESMF_GridGetCoord(make_grid, localDE=0, staggerLoc=ESMF_STAGGERLOC_CENTER, & ! coordDim=1, farrayPtr=coordX, rc=localrc) !if (ESMF_LogFoundError(localrc, & ! ESMF_ERR_PASSTHRU, & ! ESMF_CONTEXT, rcToReturn=rc)) return !do i = lbound(coordX,1), ubound(coordX,1) ! coordX(i) = startx + atm_dx/2. + (i-1)*atm_dx !enddo !print *, 'coordX: ', coordX ! X corner call ESMF_GridGetCoord(make_grid, localDE=0, staggerLoc=ESMF_STAGGERLOC_CORNER, & coordDim=1, farrayPtr=coordX, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return do i = lbound(coordX,1), ubound(coordX,1) coordX(i) = startx + (i-1)*atm_dx enddo !print *, 'startx: ', startx, lbound(coordX, 1), 'coordX: ', coordX ! Y center !call ESMF_GridGetCoord(make_grid, localDE=0, staggerLoc=ESMF_STAGGERLOC_CENTER, & ! coordDim=2, farrayPtr=coordY, rc=localrc) !if (ESMF_LogFoundError(localrc, & ! ESMF_ERR_PASSTHRU, & ! ESMF_CONTEXT, rcToReturn=rc)) return !do i = lbound(coordY,1), ubound(coordY,1) ! coordY(i) = starty + atm_dy/2. + (i-1)*atm_dy !enddo ! Y corner call ESMF_GridGetCoord(make_grid, localDE=0, staggerLoc=ESMF_STAGGERLOC_CORNER, & coordDim=2, farrayPtr=coordY, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return do i = lbound(coordY,1), ubound(coordY,1) coordY(i) = starty + (i-1)*atm_dy enddo if(present(field)) then field = ESMF_FieldCreate(make_grid, typekind=ESMF_TYPEKIND_R8, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return endif if(present(rc)) rc = ESMF_SUCCESS end function make_grid function make_grid_sph(atm_nx, atm_ny, atm_dx, atm_dy, atm_sx, atm_sy, area_adj, tag, scheme, rc) ! return value type(ESMF_Grid) :: make_grid_sph ! arguments integer, intent(in) :: atm_nx, atm_ny real(ESMF_KIND_R4), intent(in) :: atm_dx, atm_dy real(ESMF_KIND_R4), intent(in) :: atm_sx, atm_sy real(ESMF_KIND_R4), intent(in), optional :: area_adj character(len=*), intent(in), optional :: tag integer, intent(in) , optional :: scheme integer, intent(out), optional :: rc ! local variables integer :: localrc, i, j real(ESMF_KIND_R8), pointer :: coordX(:,:), coordY(:,:) real(ESMF_KIND_R8), pointer :: f_area(:,:), f_area_m(:), o_area(:,:) real(ESMF_KIND_R8) :: startx, starty integer :: l_scheme type(ESMF_Mesh) :: mesh type(ESMF_Field) :: field l_scheme = ESMF_REGRID_SCHEME_REGION3D if(present(scheme)) l_scheme = scheme if(l_scheme == ESMF_REGRID_SCHEME_FULL3D) then make_grid_sph = ESMF_GridCreate1PeriDim(maxIndex=(/atm_nx, atm_ny/), & indexflag=ESMF_INDEX_GLOBAL, & gridEdgeLWidth=(/0,0/), gridEdgeUWidth=(/0,1/), & !regDecomp=(/npet, 1/), & rc=localrc) else make_grid_sph = ESMF_GridCreateNoPeriDim(maxIndex=(/atm_nx, atm_ny/), & indexflag=ESMF_INDEX_GLOBAL, & gridEdgeLWidth=(/0,0/), gridEdgeUWidth=(/1,1/), & !regDecomp=(/npet, 1/), & rc=localrc) endif if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_GridAddCoord(make_grid_sph, staggerloc=ESMF_STAGGERLOC_CENTER, & rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_GridAddCoord(make_grid_sph, staggerloc=ESMF_STAGGERLOC_CORNER, & rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! global indexing ! atm grid is not decomposed in the y direction !startx = lpet*atm_nx/npet*atm_dx startx = atm_sx starty = atm_sy ! compute coord ! X center call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CENTER, & coordDim=1, farrayPtr=coordX, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Y center call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CENTER, & coordDim=2, farrayPtr=coordY, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return do i = lbound(coordX,1), ubound(coordX,1) do j = lbound(coordX, 2), ubound(coordX, 2) coordX(i,j) = startx + atm_dx/2. + (i-1)*atm_dx coordY(i,j) = starty + atm_dy/2. + (j-1)*atm_dy enddo enddo !print *, 'startx: ', startx, lbound(coordX, 1), ubound(coordX, 1), 'coordX: ', coordX(:,1) ! X corner call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CORNER, & coordDim=1, farrayPtr=coordX, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Y corner call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CORNER, & coordDim=2, farrayPtr=coordY, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return do i = lbound(coordX,1), ubound(coordX,1) do j = lbound(coordX, 2), ubound(coordX, 2) coordX(i,j) = startx + (i-1)*atm_dx coordY(i,j) = starty + (j-1)*atm_dy enddo enddo if(present(area_adj)) then ! retrieve area !mesh = ESMF_GridToMesh(make_grid_sph, & ! ESMF_STAGGERLOC_CORNER, 0, & ! regridConserve=ESMF_REGRID_CONSERVE_ON, rc=localrc) !if (ESMF_LogFoundError(localrc, & ! ESMF_ERR_PASSTHRU, & ! ESMF_CONTEXT, rcToReturn=rc)) return !allocate(f_area_m(mesh%NumOwnedElements)) !call ESMF_MeshGetElemArea(mesh, arealist=f_area_m, rc=localrc) !if (ESMF_LogFoundError(localrc, & ! ESMF_ERR_PASSTHRU, & ! ESMF_CONTEXT, rcToReturn=rc)) return !deallocate(f_area_m) ! find out original Grid cell area field = ESMF_FieldCreate(make_grid_sph, typekind=ESMF_TYPEKIND_R8, & staggerloc=ESMF_STAGGERLOC_CENTER, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_FieldRegridGetArea(field, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_FieldGet(field, farrayPtr=o_area, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! add area to Grid call ESMF_GridAddItem(make_grid_sph, ESMF_GRIDITEM_AREA, & staggerloc=ESMF_STAGGERLOC_CENTER, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_GridGetItem(make_grid_sph, ESMF_GRIDITEM_AREA, & staggerloc=ESMF_STAGGERLOC_CENTER, farrayptr=f_area, & rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! adjust Grid area f_area = area_adj*o_area endif if(present(rc)) rc = ESMF_SUCCESS end function make_grid_sph function make_mesh_sph(atm_nx, atm_ny, atm_dx, atm_dy, atm_sx, atm_sy, tag, scheme, rc) ! return value type(ESMF_Mesh) :: make_mesh_sph ! arguments integer, intent(in) :: atm_nx, atm_ny real(ESMF_KIND_R4), intent(in) :: atm_dx, atm_dy real(ESMF_KIND_R4), intent(in) :: atm_sx, atm_sy character(len=*), intent(in), optional :: tag integer, intent(in) , optional :: scheme integer, intent(out), optional :: rc ! local variables integer :: localrc, i, j real(ESMF_KIND_R8), pointer :: coordX(:,:), coordY(:,:) real(ESMF_KIND_R8) :: startx, starty integer :: l_scheme type(ESMF_Grid) :: make_grid_sph l_scheme = ESMF_REGRID_SCHEME_REGION3D if(present(scheme)) l_scheme = scheme if(l_scheme == ESMF_REGRID_SCHEME_FULL3D) then make_grid_sph = ESMF_GridCreate1PeriDim(maxIndex=(/atm_nx, atm_ny/), & indexflag=ESMF_INDEX_GLOBAL, & gridEdgeLWidth=(/0,0/), gridEdgeUWidth=(/0,1/), & !regDecomp=(/npet, 1/), & rc=localrc) else make_grid_sph = ESMF_GridCreateNoPeriDim(maxIndex=(/atm_nx, atm_ny/), & indexflag=ESMF_INDEX_GLOBAL, & gridEdgeLWidth=(/0,0/), gridEdgeUWidth=(/1,1/), & !regDecomp=(/npet, 1/), & rc=localrc) endif if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_GridAddCoord(make_grid_sph, staggerloc=ESMF_STAGGERLOC_CENTER, & rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_GridAddCoord(make_grid_sph, staggerloc=ESMF_STAGGERLOC_CORNER, & rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! global indexing ! atm grid is not decomposed in the y direction !startx = lpet*atm_nx/npet*atm_dx startx = atm_sx starty = atm_sy ! compute coord ! X center call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CENTER, & coordDim=1, farrayPtr=coordX, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Y center call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CENTER, & coordDim=2, farrayPtr=coordY, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return do i = lbound(coordX,1), ubound(coordX,1) do j = lbound(coordX, 2), ubound(coordX, 2) coordX(i,j) = startx + atm_dx/2. + (i-1)*atm_dx coordY(i,j) = starty + atm_dy/2. + (j-1)*atm_dy enddo enddo !print *, 'startx: ', startx, lbound(coordX, 1), ubound(coordX, 1), 'coordX: ', coordX(:,1) ! X corner call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CORNER, & coordDim=1, farrayPtr=coordX, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Y corner call ESMF_GridGetCoord(make_grid_sph, localDE=0, staggerLoc=ESMF_STAGGERLOC_CORNER, & coordDim=2, farrayPtr=coordY, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return do i = lbound(coordX,1), ubound(coordX,1) do j = lbound(coordX, 2), ubound(coordX, 2) coordX(i,j) = startx + (i-1)*atm_dx coordY(i,j) = starty + (j-1)*atm_dy enddo enddo make_mesh_sph = ESMF_GridToMesh(make_grid_sph, & ESMF_STAGGERLOC_CORNER, 0, & regridConserve=ESMF_REGRID_CONSERVE_ON, rc=localrc) if (ESMF_LogFoundError(localrc, & ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return if(present(rc)) rc = ESMF_SUCCESS end function make_mesh_sph end program ESMF_MeshOpUTest