pph_dyn parametrization added to ppmix loop

obl_inmom.f90

@@ -8,8 +8,12 @@ module obl_inmom
   subroutine init_vermix(richnum_mode, kh_km_mode, uu, vv, lu, dx, dy, dxh, dyh, hhu, hhv, hhq, zw, g, rh0, &
     rit, den, rlh, taux, tauy, aice0, tt, ss, anzt, anzu, anumaxt, anumaxu, anubgrt, anubgru)
     use obl_legacy
-    use obl_pph
-    use obl_pph_dyn
+    use obl_pph, only: pph_kh => get_eddy_diffusivity_vec, &
+                        pph_km => get_eddy_viscosity_vec, &
+                        pphParamType
+    use obl_pph_dyn, only: pph_dyn_kh => get_eddy_diffusivity_vec, &
+                            pph_dyn_km => get_eddy_viscosity_vec, &
+                            pphDynParamType
     integer, intent(in) :: richnum_mode, kh_km_mode
     real, intent(inout) :: anzt(:,:,:) ! kh
     real, intent(inout) :: anzu(:,:,:) ! km
@@ -31,6 +35,11 @@ module obl_inmom
     real, intent(in) :: anubgrt ! kh_b0
     real, intent(in) :: anubgru ! km_b0
     ! Local variables
+    integer :: nx, ny, nz
+    integer :: i, j
+
+    ! Modules types and vars
+    !! obl_legacy vars
     real, allocatable :: u2(:,:,:), v2(:,:,:)
     real, allocatable :: s2(:,:,:)
     real, allocatable :: n2(:,:,:)
@@ -43,12 +52,14 @@ module obl_inmom
     real :: kh_unstable
     real :: kh_b0, km_b0
     real :: kh_0, km_0
-    integer :: nx, ny, nz
-    integer :: i, j
-    ! Modules types
+    !! pph vars
     type(pphParamType) :: pphParams
-    type(gridDataType) :: inmomGrid
+    real :: neutral_mld(size(rit, 1), size(rit, 2))
+    real :: u_dynH(size(rit, 1), size(rit, 2))
+    !! pph_dyn vars
+    type(pphDynParamType) :: pphDynParams
 
+    !! allocate obl_legacy vars
     allocate(u2(size(uu, 1), size(uu, 2), size(uu, 3)))
     allocate(v2(size(vv, 1), size(vv, 2), size(vv, 3)))
     allocate(s2(size(u2, 1), size(u2, 2), size(u2, 3) - 1))
@@ -82,7 +93,7 @@ module obl_inmom
       call legacy_rit_top(rlh, taux, tauy, border_shift, lu, rit(:,:,1))
     end if
     ! print *, "hello aft rit"
-    call sync_xy_border(rit)
+    call sync_xy_border_3d(rit)
 
     do j = 1, ny
       do i = 1, nx
@@ -101,13 +112,16 @@ module obl_inmom
           else if (kh_km_mode == 2) then
             ! call get_Kh(kh(i,j,:), rit(i,j,:), nz)
             ! call get_Km(km(i,j,:), rit(i,j,:), nz)
-            call get_eddy_diffusivity(kh(i,j,:), rit(i,j,:), pphParams, inmomGrid)
-            call get_eddy_viscosity(km(i,j,:), rit(i,j,:), pphParams, inmomGrid)
+            call pph_kh(kh(i,j,:), rit(i,j,:), pphParams, nz)
+            call pph_km(km(i,j,:), rit(i,j,:), pphParams, nz)
           else if (kh_km_mode == 3) then
-            ! u_dynH = 
-            ! mld = 
-            call get_eddy_diffusivity(kh(i,j,:), rit(i,j,:), u_dynH, mld, pphParams, inmomGrid)
-            call get_eddy_viscosity(km(i,j,:), rit(i,j,:), u_dynH, mld, pphParams, inmomGrid)
+            if (i == 1 .and. j == 1) then
+                call legacy_neutral_mld(rlh, taux, tauy, border_shift, lu, neutral_mld, u_dynH)
+                call sync_xy_border_2d(neutral_mld)
+                call sync_xy_border_2d(u_dynH)
+            end if
+            call pph_dyn_kh(kh(i,j,:), rit(i,j,:), u_dynH(i,j), neutral_mld(i,j), pphdynParams, nz)
+            call pph_dyn_km(km(i,j,:), rit(i,j,:), u_dynH(i,j), neutral_mld(i,j), pphdynParams, nz)
           end if
         end if
       end do
@@ -127,7 +141,7 @@ module obl_inmom
     anzu = km
   end subroutine init_vermix
 
-  subroutine sync_xy_border(array)
+  subroutine sync_xy_border_3d(array)
     implicit none
     real, intent(inout) :: array(:,:,:)
     integer :: nx, ny, nz, k
@@ -157,5 +171,30 @@ module obl_inmom
       array(nx,ny,k) = (array(nx-1,ny,k) + array(nx,ny-1,k)) / 2.0
     end do
   
-  end subroutine sync_xy_border
+  end subroutine sync_xy_border_3d
+
+  subroutine sync_xy_border_2d(array)
+    implicit none
+    real, intent(inout) :: array(:,:)
+    integer :: nx, ny
+
+    ! Determine array dimensions
+    nx = size(array, 1)
+    ny = size(array, 2)
+
+    ! Update boundary points along x (first_x and end_x boundaries)
+    array(1,2:ny-1) = array(2,2:ny-1)         ! first_x boundary
+    array(nx,2:ny-1) = array(nx-1,2:ny-1)     ! end_x boundary
+
+    ! Update boundary points along y (first_y and end_y boundaries)
+    array(2:nx-1,1) = array(2:nx-1,2)         ! first_y boundary
+    array(2:nx-1,ny) = array(2:nx-1,ny-1)     ! end_y boundary
+
+    ! Update corner points
+    array(1,1) = (array(2,1) + array(1,2)) / 2.0
+    array(1,ny) = (array(2,ny) + array(1,ny-1)) / 2.0
+    array(nx,1) = (array(nx-1,1) + array(nx,2)) / 2.0
+    array(nx,ny) = (array(nx-1,ny) + array(nx,ny-1)) / 2.0
+
+  end subroutine sync_xy_border_2d
 end module obl_inmom
\ No newline at end of file