Fix issue with parallel derivatives

.github/workflows/test.yaml

@@ -0,0 +1,55 @@
+name: Test
+on: [push]
+
+jobs:
+  standard_tests:
+    name: Tests with ${{ matrix.config.name }}
+    timeout-minutes: 60
+    runs-on: ubuntu-latest
+    container:
+      image: oi4ai/bout3d:db-outer-${{ matrix.config.mode }}
+    strategy:
+      fail-fast: true
+      matrix:
+        config:
+          - name: debugging
+            mode: debug
+          - name: optimisations
+            mode: opt
+    steps:
+      - name: Job information
+        run: |
+          echo Build: ${{ matrix.config.name }}, ${{ matrix.config.os }}
+          echo CMake options: ${{ matrix.config.cmake_options }}
+          cat /etc/os-release
+
+      - name: Setup dependencies
+        run: |
+          sudo dnf install -y git-lfs
+
+      - name: Checkout hermes-2
+        run: |
+          cd /home/boutuser/
+          git clone ${{ github.server_url }}/${{ github.repository }} -b ${{ github.ref_name }} hermes-2
+          cd hermes-2
+          git checkout ${{ github.sha }}
+          git submodule update --init --recursive
+
+      - name: Build
+        run: |
+          cd /home/boutuser/hermes-2
+          cmake --version
+          cmake -S . -B build -Dbout++_DIR=/home/boutuser/BOUT-dev/build
+          make -C build -j 2
+
+      - name: Run
+        run: |
+          cd /home/boutuser/hermes-2
+          # Run only shortly
+          sed -e 's/nout = .*/nout = 1/' -i */BOUT.inp
+          sed -e 's/timestep = 100/timestep = 1/' -i */BOUT.inp
+          # Change solver type - beuler is quite heavy on setup, but we want
+          # to run only short
+          sed -e 's/type = beuler/type = pvode/' -i */BOUT.inp
+          build/hermes-2 -d 1-lowres
+          build/hermes-2 -d 3-medium

.gitignore

@@ -1,4 +1,13 @@
 *.o
 *.nc
 BOUT.log.*
-
+.BOUT.pid.*
+BOUT.settings
+/atomicpp.a
+/doc/hermes-2-manual.aux
+/doc/hermes-2-manual.log
+/doc/hermes-2-manual.out
+/doc/hermes-2-manual.pdf
+/hermes-2
+/core*
+/Makefile

.gitmodules

@@ -0,0 +1,3 @@
+[submodule "grids"]
+	path = grids
+	url = https://github.com/dschwoerer/hermes-grids.git

1-lowres/BOUT.inp

@@ -1,23 +1,23 @@
 
-NOUT = 100
-TIMESTEP = 10
+nout = 10
+timestep = 100
 
-TwistShift = false  # use twist-shift condition?
+twistshift = false  # use twist-shift condition?
 ballooning = false  # ballooning transformation?
 shiftinitial = false
-ShiftWithoutTwist=false
+ShiftWithoutTwist = false
 
-MYG=1
-MXG=2
+MYG = 1
+MXG = 2
 
-grid = "rotating-ellipse-36x8x32-ix20.fci.nc"
+grid = "grids/rotating-ellipse-36x8x32-ix20.fci.nc"
 
-NXPE=1
+NXPE = 1
 
 dump_format = nc
 
-[restart]
-init_missing=true # initialize missing variables on restart?
+[restart_files]
+init_missing = true # initialize missing variables on restart?
 
 [mesh]
 symmetricGlobalX = true
@@ -30,19 +30,19 @@ type = fci
 ##################################################
 # derivative methods
 
-[ddx]
+[mesh:ddx]
 
 first = C2
 second = C2
 upwind = W3
 
-[ddy]
+[mesh:ddy]
 
 first = C2
 second = C2
 upwind = W3
 
-[ddz]
+[mesh:ddz]
 
 first = C2
 second = C2
@@ -55,10 +55,11 @@ upwind = W3
 
 # Note: If evolving neutrals, need preconditioning
 # type = cvode
+type = beuler
 # use_precon = true
 
-ATOL = 1.0e-10  # absolute tolerance
-RTOL = 1.0e-5   # relative tolerance
+atol = 1.0e-10  # absolute tolerance
+rtol = 1.0e-5   # relative tolerance
 mxstep = 1000000  # Maximum internal steps per output
 
 cvode_max_order = 2
@@ -72,7 +73,7 @@ nonuniform = true
 # Electrostatic potential solver
 
 [phiSolver]
-type   = petsc  # Needed if Boussinesq = false
+type = petsc  # Needed if Boussinesq = false
 ksptype = gmres # Linear iterative method
 pctype = lu # Preconditioner. Direct "lu" or "ilu"; iterative "jacobi", "sor"
 
@@ -115,7 +116,8 @@ nonuniform = true
 # general settings for the model
 [input]
 
-transform_from_field_aligned=False
+transform_from_field_aligned = false
+error_on_unused_options = false
 
 [Hermes]
 loadmetric = false
@@ -140,7 +142,7 @@ resistivity_boundary_width = 0
 # Numerical dissipation
 
 vepsi_dissipation = false  # Parallel dissipation on Ve-Vi
-vort_dissipation = false
+vort_dissipation = true # false # Changed 2021-10-25
 numdiff = 0.01
 hyperpar = -1
 
@@ -171,16 +173,18 @@ split_n0_psi = false
 j_diamag       = false  # Diamagnetic current: Vort <-> Pe
 j_par          = false  # Parallel current:    Vort <-> Psi
 
-evolve_te = false
+evolve_te = true
 evolve_ti = false
 
+evolve_vort = false
+
 pe_par         = true  # Parallel pressure gradient: Pe <-> Psi
 resistivity    = true  # Resistivity: Psi -> Pe
 thermal_flux   = false
 thermal_force  = false
 electron_ion_transfer = false
 electron_viscosity = false
-ion_viscosity  = false  # Ion parallel viscosity
+ion_viscosity = false  # Ion parallel viscosity
 thermal_conduction = false
 
 frecycle = 0.0  # Neutral gas recycling fraction
@@ -199,16 +203,16 @@ drift_wave = false
 # Transport coefficients
 classical_diffusion = false  # Collisional diffusion
 
-anomalous_D = -1  # Anomalous density diffusion [m^2/s]
-anomalous_chi = -1 # Anomalous thermal diffusion [m^2/s]
+anomalous_D = 0.2  # Anomalous density diffusion [m^2/s]
+anomalous_chi = 0.6 # Anomalous thermal diffusion [m^2/s]
 anomalous_nu = -1   # Anomalous viscosity
 
 poloidal_flows = false
 
 magnetic_drift = true
 ion_velocity   = true
 
-ion_neutral = 0.0
+ion_neutral = false
 neutral_friction = false  # Friction between plasma and neutrals
 
 boussinesq = true   # Use Boussinesq approximation
@@ -229,6 +233,7 @@ sheath_ydown = true
 sheath_allow_supersonic = false
 sheath_gamma_e = 4   # Electron sheath heat transmission
 sheath_gamma_i = 2.5 # Ion sheath heat transmission
+parallel_sheaths = true
 neutral_gamma = 0.0
 
 startprofiles = false
@@ -254,18 +259,21 @@ AA = 1  # Atomic mass. 1 = Hydrogen, 2 = Deuterium
 type = none    # Neutral model: none, diffusion2d, recycling, fullvelocity, mixed
 viscosity = 1  # Dynamic viscosity
 bulk = 0       # Bulk (volume) viscosity
-conduction = 1 
+conduction = 1
 neutral_gamma = 0.0
 
 nn_floor = 1e-2  # Floor applied when calculating Vn = NVn / Nn
 low_n_equilibriate = -1e-4  # If n < nn_floor, evolve Tn and Vn towards plasma values
 
-[All]
+[all]
 scale = 0.0
 
 bndry_all = neumann_o2
 bndry_xin = neumann_o2
 bndry_xout = neumann_o2
+# Maybe not documented?
+bndry_par_all = parallel_neumann
+
 
 [Ne] # Electron density
 scale = 1
@@ -275,9 +283,11 @@ function = 0.1 + 0.01*gauss(x,0.21)
 #function = 1 + 0.1*gauss(x-0.5, 0.14)*gauss(z-pi/2.,0.14/2.) + 0.1*gauss(x-0.5, 0.14)*gauss(z-pi,0.14/2.) + 0.1*gauss(x-0.5, 0.14)*gauss(z,0.14/2.) + 0.1*gauss(x-0.5, 0.14)*gauss(z-3*pi/2.,0.14/2.)
 
 source = 2.8e3*gauss(x, 0.21)
+bndry_xin = neumann_o2
+#dirichlet(0.1)
 
 [Vort]
-scale=0
+scale = 0
 function = sin(7*z)*sin(3*x-z)*gauss(-x,0.21)*(1+sin(y))#mixmode(y)
 bndry_all = dirichlet_o2
 
@@ -291,20 +301,22 @@ scale = 1
 function = 0.1 + 0.01*cos(z)*(sin(8*pi*x - 30*z)+sin(7*z)+0.5*sin(13*z - 2*pi*x))*gauss(-x,0.21)
 #function = 1 + 1.2*(Ne:function-1)
 source = Ne:source
+# bndry_par_all = parallel_free
 
 [Pi]
 scale = 1
 #function = 1 + 1.2*(Ne:function-1)
 function = 0.1 + 0.01*cos(z)*(sin(8*pi*x - 30*z)+sin(7*z)+0.5*sin(13*z - 2*pi*x))*gauss(-x,0.21)
 source = Ne:source
+# bndry_xin = dirichlet(1)
 
 [Ve]
-parallel_bndry_yup=1
-parallel_bndry_ydown=-1
+parallel_bndry_yup = 1
+parallel_bndry_ydown = -1
 
 [Jpar]
-parallel_bndry_yup=1
-parallel_bndry_ydown=-1
+parallel_bndry_yup = 1
+parallel_bndry_ydown = -1
 
 [phi]
 # Radial boundaries determined by Laplacian inversion
@@ -348,4 +360,3 @@ bndry_all = neumann_o2
 
 [NVn]
 bndry_all = dirichlet_o2
-