From 4b17b50d6e0c3c7bf6d8bb72fd25af187dfb522e Mon Sep 17 00:00:00 2001
From: gianlupo <40678236+gianlupo@users.noreply.github.com>
Date: Thu, 18 Jul 2024 18:48:47 +0200
Subject: [PATCH] 3D Antuono vortex initial field. (#123)

---
 docs/INFO_INPUT.md |  1 +
 src/initflow.f90   | 25 ++++++++++++++++++++++++-
 2 files changed, 25 insertions(+), 1 deletion(-)

diff --git a/docs/INFO_INPUT.md b/docs/INFO_INPUT.md
index 39f6773d..9768ef82 100644
--- a/docs/INFO_INPUT.md
+++ b/docs/INFO_INPUT.md
@@ -106,6 +106,7 @@ These lines set the initial velocity field.
 * `hdc`: half plane Poiseuille flow profile with constant pressure gradient         ; streamwise direction in `x`
 * `tgv`: three-dimensional Taylor-Green vortex
 * `tgw`: two-dimensional   Taylor-Green vortex
+* `ant`: three-dimensional Antuono vortex
 
 `is_wallturb`, if true, **superimposes a high amplitude disturbance on the initial velocity field** that effectively triggers transition to turbulence in a wall-bounded shear flow.
 
diff --git a/src/initflow.f90 b/src/initflow.f90
index 0f36a2dd..5e3cc24e 100644
--- a/src/initflow.f90
+++ b/src/initflow.f90
@@ -114,6 +114,29 @@ subroutine initflow(inivel,bcvel,ng,lo,l,dl,zc,zf,dzc,dzf,visc, &
           end do
         end do
       end do
+    case('ant')
+      !
+      ! see M. Antuono, JFM 890, A23 (2020)
+      !
+      do k=1,n(3)
+        zcc = zc(k)/l(3)*2.*pi+0.5*pi
+        zff = zf(k)/l(3)*2.*pi+0.5*pi
+        do j=1,n(2)
+          yc = (j+lo(2)-1-.5)*dl(2)/l(2)*2.*pi+0.5*pi
+          yf = (j+lo(2)-1-.0)*dl(2)/l(2)*2.*pi+0.5*pi
+          do i=1,n(1)
+            xc = (i+lo(1)-1-.5)*dl(1)/l(1)*2.*pi+0.5*pi
+            xf = (i+lo(1)-1-.0)*dl(1)/l(1)*2.*pi+0.5*pi
+            u(i,j,k) = (4.*sqrt(2.)/3./sqrt(3.))*(sin(xf-5.*pi/6.)*cos(yc-1.*pi/6.)*sin(zcc         ) - &
+                                                  sin(xf-1.*pi/6.)*sin(yc         )*cos(zcc-5.*pi/6.))*uref
+            v(i,j,k) = (4.*sqrt(2.)/3./sqrt(3.))*(sin(xc         )*sin(yf-5.*pi/6.)*sin(zcc-1.*pi/6.) - &
+                                                  cos(xc-5.*pi/6.)*sin(yf-1.*pi/6.)*sin(zcc         ))*uref
+            w(i,j,k) = (4.*sqrt(2.)/3./sqrt(3.))*(cos(xc-1.*pi/6.)*sin(yc         )*sin(zff-5.*pi/6.) - &
+                                                  sin(xc         )*cos(yc-5.*pi/6.)*sin(zff-1.*pi/6.))*uref
+            p(i,j,k) = -(u(i,j,k)**2+v(i,j,k)**2+w(i,j,k)**2)/2.
+          end do
+        end do
+      end do
     case('tgw')
       do k=1,n(3)
         do j=1,n(2)
@@ -161,7 +184,7 @@ subroutine initflow(inivel,bcvel,ng,lo,l,dl,zc,zf,dzc,dzf,visc, &
       call MPI_FINALIZE(ierr)
       error stop
     end select
-    if(.not.any(inivel == ['tgv','tgw'])) then
+    if(.not.any(inivel == ['tgv','tgw','ant'])) then
       do k=1,n(3)
         do j=1,n(2)
           do i=1,n(1)