Rotordisk class and stirred tank example

examples/stirredTank/GNUmakefile

@@ -0,0 +1,50 @@
+# NGA location if not yet defined
+NGA_HOME ?= ../..
+
+# Compilation parameters
+PRECISION = DOUBLE
+USE_MPI   = TRUE
+USE_HYPRE = TRUE
+USE_LAPACK= TRUE
+USE_IRL   = FALSE
+USE_FFTW  = TRUE
+PROFILE   = FALSE
+DEBUG     = FALSE
+COMP      = gnu
+EXEBASE   = nga
+
+# Directories that contain user-defined code
+Udirs   := src
+
+# Include user-defined sources
+Upack   += $(foreach dir, $(Udirs), $(wildcard $(dir)/Make.package))
+Ulocs   += $(foreach dir, $(Udirs), $(wildcard $(dir)))
+include $(Upack)
+INCLUDE_LOCATIONS += $(Ulocs)
+VPATH_LOCATIONS   += $(Ulocs)
+
+# External libraries are defined in .profile/.bashrc/.zshrc, but could be defined here as well
+
+# NGA compilation definitions
+include $(NGA_HOME)/tools/GNUMake/Make.defs
+
+# Include NGA base code
+Bdirs   := variable_density core data solver config grid libraries rotor_disk transform subgrid 
+Bpack   += $(foreach dir, $(Bdirs), $(NGA_HOME)/src/$(dir)/Make.package)
+include $(Bpack)
+
+# Inform user of Make.packages used
+ifdef Ulocs
+   $(info Taking user code from: $(Ulocs))
+endif
+$(info Taking base code from: $(Bdirs))
+
+# Target definition
+all: $(executable)
+	@echo COMPILATION SUCCESSFUL
+
+# NGA compilation rules
+include $(NGA_HOME)/tools/GNUMake/Make.rules
+
+run:
+	mpiexec -n 8 $(executable) -i input -v 2

examples/stirredTank/input

@@ -0,0 +1,45 @@
+# Parallelization
+Partition : 1 8 1
+
+# Mesh definition
+Lx : 0.32
+Ly : 0.32
+Lz : 0.32
+nx : 96
+ny : 96
+nz : 96
+tank radius : 0.14
+tank height : 0.28
+
+# Fluid properties
+Gravity :                   0 -9.81 0
+Dynamic viscosity:          0.001
+Density :                   1000
+
+# Rotor Disk properties
+# Blade
+Number of radical points : 2
+Blade radius : 0.03 0.08
+Blade twists : -30.0 -30.0
+Blade Chords : 0.02 0.02
+Number of azimuthal points : 9
+# estimation of plate airfoil
+Blade AoA : -180 -135 -90 -45 0 45 90 135 180
+Blade Cd : 0 1.0 1.25 1.0 0 1.0 1.25 1.0 0
+Blade Cl : 0 1.0 0 -1.0 0 1.0 0 -1.0 0
+# Rotor Disk
+Rotor Disk number of blades : 4
+Rotor Disk min radius : 0.03    # matching the blade radius
+Rotor Disk max radius : 0.08
+Rotor Disk center : 0.0 0.08 0.0
+Rotor Disk axis : 0.0 1.0 0.0
+Rotor Disk reference direction : 1.0 0.0 0.0  # perpendicular to the axis direction
+Rotor Disk angular velocity : 120      # rpm
+
+# Time integration
+Max timestep size : 2e-3
+Max cfl number :    0.9
+Max time :          5
+
+# Ensight output
+Ensight output period : 5e-2

examples/stirredTank/src/Make.package

@@ -0,0 +1,2 @@
+# List here the extra files here
+f90EXE_sources += geometry.f90 simulation.f90

examples/stirredTank/src/geometry.f90

@@ -0,0 +1,107 @@
+!> Various definitions and tools for initializing NGA2 config
+module geometry
+   use ibconfig_class, only: ibconfig
+   use precision,    only: WP
+   implicit none
+   private
+
+   !> Single config
+   type(ibconfig), public :: cfg
+
+   !> Tank radius
+   real(WP), public :: r_tank
+   real(WP), public :: height_tank
+
+   public :: geometry_init
+
+contains
+
+
+   !> Initialization of problem geometry
+   subroutine geometry_init
+      use sgrid_class, only: sgrid
+      use param,       only: param_read
+      implicit none
+      type(sgrid) :: grid
+
+      ! Create a grid from input params
+      create_grid: block
+         use sgrid_class, only: cartesian
+         integer :: i,j,k,nx,ny,nz
+         real(WP) :: Lx,Ly,Lz
+         real(WP), dimension(:), allocatable :: x,y,z
+
+         ! Read in grid definition
+         call param_read('Lx',Lx); call param_read('nx',nx); allocate(x(nx+1))
+         call param_read('Ly',Ly); call param_read('ny',ny); allocate(y(ny+1))
+         call param_read('Lz',Lz); call param_read('nz',nz); allocate(z(nz+1))
+
+         ! Create simple rectilinear grid
+         do i=1,nx+1
+            x(i)=real(i-1,WP)*Lx/real(nx,WP)-0.5_WP*Lx
+         end do
+         do j=1,ny+1
+            y(j)=real(j-1,WP)*Ly/real(ny,WP)
+         end do
+         do k=1,nz+1
+            z(k)=real(k-1,WP)*Lz/real(nz,WP)-0.5_WP*Lz
+         end do
+
+         ! General serial grid object
+         grid=sgrid(coord=cartesian,no=3,x=x,y=y,z=z,xper=.true.,yper=.true.,zper=.true.,name='Stirred Tank')
+
+      end block create_grid
+
+      ! Create a config from that grid on our entire group
+      create_cfg: block
+         use parallel, only: group
+         integer, dimension(3) :: partition
+
+         ! Read in partition
+         call param_read('Partition',partition,short='p')
+
+         ! Create partitioned grid
+         cfg=ibconfig(grp=group,decomp=partition,grid=grid)
+
+      end block create_cfg
+
+     ! Create masks for this config
+      create_walls: block
+         use ibconfig_class, only: bigot,sharp
+         use param,          only: param_read
+         integer :: i,j,k
+         real(WP) :: xm, ym, zm, tmp1, tmp2
+         real(WP) :: Ly
+
+         ! Read in tank radius
+         call param_read('tank radius',r_tank)
+         call param_read('tank height',height_tank)
+         call param_read('Ly',Ly)
+
+         do k=cfg%kmino_,cfg%kmaxo_
+            do j=cfg%jmino_,cfg%jmaxo_
+               do i=cfg%imino_,cfg%imaxo_
+                  xm=cfg%xm(i); ym=cfg%ym(j); zm=cfg%zm(k)
+                  tmp1 = sqrt(xm**2 + zm**2) - r_tank ! distance to wall
+                  if (ym .gt. 0.5_WP*(Ly+height_tank)) then
+                     tmp2 = ym - 0.5_WP*(Ly+height_tank) ! distance to top
+                  else if (ym .lt. 0.5_WP*(Ly-height_tank)) then
+                     tmp2 = 0.5_WP*(Ly-height_tank) - ym
+                  else
+                     tmp2 = min(ym-0.5_WP*(Ly+height_tank),0.5_WP*(Ly-height_tank)-ym)
+                  end if
+                  cfg%Gib(i,j,k)=max(tmp1,tmp2)
+               end do
+            end do
+         end do
+          ! Get normal vector
+         call cfg%calculate_normal()
+         ! Get VF field
+         call cfg%calculate_vf(method=sharp,allow_zero_vf=.false.)
+      end block create_walls
+
+
+   end subroutine geometry_init
+
+
+end module geometry