[WIP] Random waves vegetation

2d/waves_vegetation/README.rst

@@ -0,0 +1,21 @@
+Flow Through Vegetation (Anderson and Smith, 2014)
+==================================================
+
+The problem's domain is a 63.4m long, 1.5m wide flume of varying depth. A
+wave is generated at the far edge of the flume, which is at a depth of
+1.95m. 
+The deep section is 5.4m-long, followed by a 1:44 slope for 19.5m
+that connects to a 12.2-m long, 1.5m-deep flat testing area. 
+Within the testing area there are rods simulating idealized vegetation placed in a diamond formation. 
+After the testing area there is a 11.3m, 1:20 slope with damping properties (to reduce reflected waves). 
+
+.. figure:: ./image_here.bmp
+   :width: 100%
+   :align: center
+
+This study uses PROTEUS to model wave attenuation in shallow-water vegetated areas. 
+
+References
+----------
+
+- Anderson, M.E. and Smith, J.M. (2013), Wave attenuation by flexible, idealized salt marsh vegetation. Coastal Engineering, Volume 83 p.82-92.

2d/waves_vegetation/lonestar-runs/ls_consrv_n.py

@@ -0,0 +1,48 @@
+from proteus import *
+from tank import *
+from ls_consrv_p import *
+
+timeIntegrator  = ForwardIntegrator
+timeIntegration = NoIntegration
+
+femSpaces = {0:basis}
+
+subgridError      = None
+massLumping       = False
+numericalFluxType = DoNothing
+conservativeFlux  = None
+shockCapturing    = None
+
+fullNewtonFlag = True
+multilevelNonlinearSolver = Newton
+levelNonlinearSolver      = Newton
+
+nonlinearSmoother = None
+linearSmoother    = None
+
+matrix = SparseMatrix
+
+if useOldPETSc:
+    multilevelLinearSolver = PETSc
+    levelLinearSolver      = PETSc
+else:
+    multilevelLinearSolver = KSP_petsc4py
+    levelLinearSolver      = KSP_petsc4py
+
+if useSuperlu:
+    multilevelLinearSolver = LU
+    levelLinearSolver      = LU
+
+linear_solver_options_prefix = 'mcorr_'
+nonlinearSolverConvergenceTest = 'r'
+levelNonlinearSolverConvergenceTest = 'r'
+linearSolverConvergenceTest  = 'r-true'
+
+tolFac = 0.0
+linTolFac = 0.01
+l_atol_res = 0.01*mcorr_nl_atol_res
+nl_atol_res = mcorr_nl_atol_res
+useEisenstatWalker = False
+
+maxNonlinearIts = 50
+maxLineSearches = 0

2d/waves_vegetation/lonestar-runs/ls_consrv_p.py

@@ -0,0 +1,26 @@
+from proteus import *
+from proteus.default_p import *
+from tank import *
+from proteus.mprans import MCorr
+
+LevelModelType = MCorr.LevelModel
+
+coefficients = MCorr.Coefficients(LSModel_index=2,V_model=0,me_model=4,VOFModel_index=1,
+                                  applyCorrection=applyCorrection,nd=nd,checkMass=False,useMetrics=useMetrics,
+                                  epsFactHeaviside=epsFact_consrv_heaviside,
+                                  epsFactDirac=epsFact_consrv_dirac,
+                                  epsFactDiffusion=epsFact_consrv_diffusion)
+
+class zero_phi:
+    def __init__(self):
+        pass
+    def uOfX(self,X):
+        return 0.0
+    def uOfXT(self,X,t):
+        return 0.0
+
+initialConditions  = {0:zero_phi()}
+
+
+
+

2d/waves_vegetation/lonestar-runs/ls_n.py

@@ -0,0 +1,62 @@
+from proteus import *
+from ls_p import *
+
+if timeDiscretization=='vbdf':
+    timeIntegration = VBDF
+    timeOrder=2
+    stepController  = Min_dt_cfl_controller
+elif timeDiscretization=='flcbdf':
+    timeIntegration = FLCBDF
+    #stepController = FLCBDF_controller
+    stepController  = Min_dt_cfl_controller
+    time_tol = 10.0*ls_nl_atol_res
+    atol_u = {0:time_tol}
+    rtol_u = {0:time_tol}
+else:
+    timeIntegration = BackwardEuler_cfl
+    stepController  = Min_dt_cfl_controller
+
+femSpaces = {0:basis}
+
+massLumping       = False
+conservativeFlux  = None
+numericalFluxType = NCLS.NumericalFlux
+subgridError      = NCLS.SubgridError(coefficients,nd)
+shockCapturing    = NCLS.ShockCapturing(coefficients,nd,shockCapturingFactor=ls_shockCapturingFactor,lag=ls_lag_shockCapturing)
+
+fullNewtonFlag  = True
+multilevelNonlinearSolver = Newton
+levelNonlinearSolver      = Newton
+
+nonlinearSmoother = None
+linearSmoother    = None
+
+matrix = SparseMatrix
+
+if useOldPETSc:
+    multilevelLinearSolver = PETSc
+    levelLinearSolver      = PETSc
+else:
+    multilevelLinearSolver = KSP_petsc4py
+    levelLinearSolver      = KSP_petsc4py
+
+if useSuperlu:
+    multilevelLinearSolver = LU
+    levelLinearSolver      = LU
+
+linear_solver_options_prefix = 'ncls_'
+nonlinearSolverConvergenceTest = 'r'
+levelNonlinearSolverConvergenceTest = 'r'
+linearSolverConvergenceTest         = 'r-true'
+
+tolFac = 0.0
+nl_atol_res = ls_nl_atol_res
+
+linTolFac = 0.0
+l_atol_res = 0.1*ls_nl_atol_res
+
+useEisenstatWalker = False
+
+maxNonlinearIts = 50
+maxLineSearches = 0
+

2d/waves_vegetation/lonestar-runs/ls_p.py

@@ -0,0 +1,29 @@
+from proteus import *
+from proteus.default_p import *
+from tank import *
+from proteus.mprans import NCLS
+
+LevelModelType = NCLS.LevelModel
+
+coefficients = NCLS.Coefficients(V_model=0,RD_model=3,ME_model=2,
+                                 checkMass=False, useMetrics=useMetrics,
+                                 epsFact=epsFact_consrv_heaviside,sc_uref=ls_sc_uref,sc_beta=ls_sc_beta,movingDomain=movingDomain)
+
+def getDBC_ls(x,flag):
+    if flag == boundaryTags['left']:
+        return wavePhi
+#    elif flag == boundaryTags['right']:
+#        return  outflowPhi
+    else:
+        return None
+
+dirichletConditions = {0:getDBC_ls}
+
+advectiveFluxBoundaryConditions =  {}
+diffusiveFluxBoundaryConditions = {0:{}}
+
+class PerturbedSurface_phi:       
+    def uOfXT(self,x,t):
+        return signedDistance(x)
+
+initialConditions  = {0:PerturbedSurface_phi()}

2d/waves_vegetation/lonestar-runs/makeRuns.py

@@ -0,0 +1,22 @@
+import numpy as np
+import os
+import shutil
+
+params = np.loadtxt("tank_parameters.csv", delimiter=',', skiprows=1)
+os.system("mkdir runs")
+for i in range(params.shape[0]):
+    filepath = "runs/"+ "run" + `i`
+    os.system("mkdir " + filepath)
+    os.system("cp *.py " + filepath)
+    os.system("cp tank.stampede.slurm " + filepath)
+    f = open(filepath + "/context.options",'w')
+    f.write("parallel=True ")
+    f.write("wave_type='single-peaked' ")
+    f.write("gauges=True ")
+    f.write("depth=" + `params[i][0]` + " ")
+    f.write("wave_height=" + `params[i][1]` + " ")
+    f.write("peak_period=" + `params[i][2]` + " ")
+    f.write("peak_wavelength=" + `params[i][3]` + " ")
+    f.write("tank_height=" + `params[i][4]`)
+    f.close()
+

2d/waves_vegetation/lonestar-runs/petsc.options.asm

@@ -0,0 +1,7 @@
+-rans2p_ksp_type gmres -rans2p_pc_type asm -rans2p_pc_asm_type basic -rans2p_ksp_max_it 2000
+-rans2p_ksp_gmres_modifiedgramschmidt -rans2p_ksp_gmres_restart 300 -rans2p_sub_ksp_type preonly -rans2p_sub_pc_factor_mat_solver_package superlu -rans2p_ksp_knoll -rans2p_sub_pc_type lu
+-ncls_ksp_type   gmres -ncls_pc_type   hypre -ncls_pc_hypre_type    boomeramg -ncls_ksp_gmres_restart 300 -ncls_ksp_knoll -ncls_ksp_max_it 2000
+-vof_ksp_type    gmres -vof_pc_type    hypre -vof_pc_hypre_type     boomeramg -vof_ksp_gmres_restart 300 -vof_ksp_knoll -vof_ksp_max_it 2000
+-rdls_ksp_type   gmres -rdls_pc_type asm -rdls_pc_asm_type basic -rdls_ksp_gmres_modifiedgramschmidt -rdls_ksp_gmres_restart 300 -rdls_ksp_knoll -rdls_sub_ksp_type preonly -rdls_sub_pc_factor_mat_solver_package superlu -rdls_sub_pc_type lu -rdls_ksp_max_it 2000
+-mcorr_ksp_type  cg -mcorr_pc_type hypre -mcorr_pc_hypre_type boomeramg -mcorr_ksp_max_it 2000
+-log_summary

2d/waves_vegetation/lonestar-runs/redist_n.py

@@ -0,0 +1,67 @@
+from proteus import *
+from redist_p import *
+from tank import *
+
+nl_atol_res = rd_nl_atol_res
+tolFac = 0.0
+nl_atol_res = rd_nl_atol_res
+
+linTolFac = 0.01
+l_atol_res = 0.01*rd_nl_atol_res
+
+if redist_Newton:
+    timeIntegration = NoIntegration
+    stepController = Newton_controller
+    maxNonlinearIts = 50
+    maxLineSearches = 0
+    nonlinearSolverConvergenceTest = 'r'
+    levelNonlinearSolverConvergenceTest = 'r'
+    linearSolverConvergenceTest = 'r-true'
+    useEisenstatWalker = False
+else:
+    timeIntegration = BackwardEuler_cfl
+    stepController = RDLS.PsiTC
+    runCFL=2.0
+    psitc['nStepsForce']=3
+    psitc['nStepsMax']=50
+    psitc['reduceRatio']=2.0
+    psitc['startRatio']=1.0
+    rtol_res[0] = 0.0
+    atol_res[0] = rd_nl_atol_res
+    useEisenstatWalker = False
+    maxNonlinearIts = 1
+    maxLineSearches = 0
+    nonlinearSolverConvergenceTest = 'rits'
+    levelNonlinearSolverConvergenceTest = 'rits'
+    linearSolverConvergenceTest = 'r-true'
+
+femSpaces = {0:basis}
+
+massLumping       = False
+numericalFluxType = DoNothing    
+conservativeFlux  = None
+subgridError      = RDLS.SubgridError(coefficients,nd)
+shockCapturing    = RDLS.ShockCapturing(coefficients,nd,shockCapturingFactor=rd_shockCapturingFactor,lag=rd_lag_shockCapturing)
+
+fullNewtonFlag = True
+multilevelNonlinearSolver  = Newton
+levelNonlinearSolver       = Newton
+
+nonlinearSmoother = NLGaussSeidel
+linearSmoother    = None
+
+matrix = SparseMatrix
+
+if useOldPETSc:
+    multilevelLinearSolver = PETSc
+    levelLinearSolver      = PETSc
+else:
+    multilevelLinearSolver = KSP_petsc4py
+    levelLinearSolver      = KSP_petsc4py
+
+if useSuperlu:
+    multilevelLinearSolver = LU
+    levelLinearSolver      = LU
+
+linear_solver_options_prefix = 'rdls_'
+

2d/waves_vegetation/lonestar-runs/redist_p.py

@@ -0,0 +1,32 @@
+from proteus import *
+from proteus.default_p import *
+from math import *
+from tank import *
+from proteus.mprans import RDLS
+"""
+The redistancing equation in the sloshbox test problem.
+"""
+
+LevelModelType = RDLS.LevelModel
+
+coefficients = RDLS.Coefficients(applyRedistancing=applyRedistancing,
+                                 epsFact=epsFact_redistance,
+                                 nModelId=2,
+                                 rdModelId=3,
+                                 useMetrics=useMetrics,
+                                 backgroundDiffusionFactor=backgroundDiffusionFactor)
+
+def getDBC_rd(x,flag):
+    pass
+
+dirichletConditions     = {0:getDBC_rd}
+weakDirichletConditions = {0:RDLS.setZeroLSweakDirichletBCsSimple}
+
+advectiveFluxBoundaryConditions =  {}
+diffusiveFluxBoundaryConditions = {0:{}}
+
+class PerturbedSurface_phi:       
+    def uOfXT(self,x,t):
+        return signedDistance(x)
+
+initialConditions  = {0:PerturbedSurface_phi()}