symm_ed Subroutine

private subroutine symm_ed(im, lamda, theta)

Arguments

Type IntentOptional Attributes Name
integer :: im
real(kind=ESMF_KIND_R8) :: lamda(im+1,im+1)
real(kind=ESMF_KIND_R8) :: theta(im+1,im+1)

Called by

proc~~symm_ed~~CalledByGraph proc~symm_ed symm_ed proc~gnomonic_grids gnomonic_grids proc~gnomonic_grids->proc~symm_ed proc~esmf_utilcreatecscoords ESMF_UtilCreateCSCoords proc~esmf_utilcreatecscoords->proc~gnomonic_grids proc~esmf_utilcreatecscoordspar ESMF_UtilCreateCSCoordsPar proc~esmf_utilcreatecscoordspar->proc~gnomonic_grids proc~esmf_gridcreatecubedsphereireg ESMF_GridCreateCubedSphereIReg proc~esmf_gridcreatecubedsphereireg->proc~esmf_utilcreatecscoordspar proc~esmf_gridcreatecubedspherereg ESMF_GridCreateCubedSphereReg proc~esmf_gridcreatecubedspherereg->proc~esmf_utilcreatecscoordspar proc~esmf_meshcreatecubedsphere ESMF_MeshCreateCubedSphere proc~esmf_meshcreatecubedsphere->proc~esmf_utilcreatecscoords interface~esmf_gridcreatecubedsphere ESMF_GridCreateCubedSphere interface~esmf_gridcreatecubedsphere->proc~esmf_gridcreatecubedsphereireg interface~esmf_gridcreatecubedsphere->proc~esmf_gridcreatecubedspherereg proc~test_bilinear_regrid_csmesh test_bilinear_regrid_csmesh proc~test_bilinear_regrid_csmesh->proc~esmf_meshcreatecubedsphere proc~test_conserve_regrid_csmesh test_conserve_regrid_csmesh proc~test_conserve_regrid_csmesh->proc~esmf_meshcreatecubedsphere proc~test_nearest_regrid_csmesh test_nearest_regrid_csmesh proc~test_nearest_regrid_csmesh->proc~esmf_meshcreatecubedsphere proc~test_patch_regrid_csmesh test_patch_regrid_csmesh proc~test_patch_regrid_csmesh->proc~esmf_meshcreatecubedsphere program~esmf_meshex ESMF_MeshEx program~esmf_meshex->proc~esmf_meshcreatecubedsphere

Source Code

 subroutine symm_ed(im, lamda, theta)
! Make grid symmetrical to i=im/2+1
 integer im
 real(ESMF_KIND_R8) lamda(im+1,im+1)
 real(ESMF_KIND_R8) theta(im+1,im+1)
 integer i,j,ip,jp
 real(ESMF_KIND_R8) avg

 do j=2,im+1
    do i=2,im
       lamda(i,j) = lamda(i,1)
    enddo
 enddo

 do j=1,im+1
    do i=1,im/2
       ip = im + 2 - i
       avg = 0.5*(lamda(i,j)-lamda(ip,j))
       lamda(i, j) = avg + pi
       lamda(ip,j) = pi - avg 
       avg = 0.5*(theta(i,j)+theta(ip,j))
       theta(i, j) = avg
       theta(ip,j) = avg
    enddo
 enddo

! Make grid symmetrical to j=im/2+1
 do j=1,im/2
       jp = im + 2 - j
    do i=2,im
       avg = 0.5*(lamda(i,j)+lamda(i,jp))
       lamda(i, j) = avg
       lamda(i,jp) = avg
       avg = 0.5*(theta(i,j)-theta(i,jp))
       theta(i, j) =  avg
       theta(i,jp) = -avg
    enddo
 enddo

 end subroutine symm_ed