Boundary Condition : PEC (#1767)

* Read boundary and set periodicity, enumerate BC types, added support for periodic

* separate particle and field boudnary structs

* PEC Etangential Bnormal and setting default pml=0

* eol

* add PEC and FieldBoundary routines

* eol

* remove duplicate struct

* update guard cells for PEC

* host device functions

* remove duplicate function

* fix compilation error by adding ;

* temporary logic to set PEC boundary

* eol

* add documentation for PEC and abort if PEC is used for RZ

* documentation for PEC functions

* remove print statements

* Adding field CI test for PEC for 3D with and without MR

* fix eol

* Apply suggestions from code review

* Update Docs/source/usage/parameters.rst

* Update Source/BoundaryConditions/WarpX_PEC.cpp

* particle PEC CI test

* rename analysis scripts for PEC

* use iside == 0 as condition

* ParallelFor over guard cells. Set field values on boundary and guard in the kernel

* ijk_guard not needed anymore

* remove unnecessary nbody generated by yt

* pi and clight are built-in

* cleanup

* eol

* fix bug in computing ig and ijk_valid

* bug fix in guard-cell update in PEC changes benchmark for particles

* fix bug in the ig compute

* reset particle PEC benchmark

* fix logic to set PEC with current interface

* adding a none type for RZ

* reset benchmarks for 2dLaserInjection, 2dgalilean_hybrid and 2dcomoving CI tests due to differences in PEC reflection

* ensure noz shape factor is used for RZ and k=0 for RZ

* fix XZ/RZ logic for Ey,Etheta,By,Btheta. Add ncomponent for RZ, generically.

* fix the field update for corner cases and mixed boundaries. rename ijk_valid to ijk_mirror since ijk_mirror need not be in valid region

* set GuardCell and Boundary cell flags

* fix compilation bug

* suppress unused variable warning

* eol whitespace

* fix doc

* Apply suggestions from code review

Co-authored-by: Andy Nonaka <AJNonaka@lbl.gov>

* RZ abort message

* reset benchmark after fixing guard-cell updates

* fix typo

* RZ doc

* Update Source/Utils/WarpXUtil.cpp

* Update Examples/Tests/PEC/inputs_particle_PEC_3d

* changes suggested from code-review

* set Enormal and Btangential in guard cell across PEC boundary

* add more info about guard cell update in the docs

* fix typo

* ijk mirror only if ig>0

* update benchmarks

* ignore unused k

* temporarily revert guardcell damping for Cartesian

* fix typo

* dont set Etangential and Bnormal to zero if field is not nodal on boundary

* set bindary int as boolean

* fix high to hi

* set default nox to 1

* particle shape for PEC test

* update benchmarks for 2D comoving and galilean

* reset benchmark for RZ test-cases

* Add in benchmarks generated by @EZoni

* Changing transverse boundaries to periodic for psatd galilean cases

* resetting benchmark for comoving_hybrid_2d and galilean_hybrid_2d

* at-least one guard cell must be filled for PEC

* abort if PEC is used for PSATD

* use new boundary interface to set none for rmin and rmax, since pec does not work for psatd

* Update Source/WarpX.cpp

* add doc mentioning PEC does not work for PSATD

* missing semicolon

* Update Source/BoundaryConditions/WarpX_PEC.H

Co-authored-by: Edoardo Zoni <59625522+EZoni@users.noreply.github.com>

* Update Source/BoundaryConditions/WarpX_PEC.H

* lengthy comment for normal and tangential components at domain boundaries

* eol fix

* grammar

* field_hi

* reset benchmarks

Co-authored-by: Andy Nonaka <AJNonaka@lbl.gov>
Co-authored-by: Edoardo Zoni <59625522+EZoni@users.noreply.github.com>

1 files changed, 84 insertions, 0 deletions

diff --git a/Examples/Tests/PEC/analysis_pec_mr.py b/Examples/Tests/PEC/analysis_pec_mr.py
new file mode 100755
index 000000000..067f2e658
--- /dev/null
+++ b/Examples/Tests/PEC/analysis_pec_mr.py
@@ -0,0 +1,84 @@
+#! /usr/bin/env python
+
+#
+#
+# This file is part of WarpX.
+#
+# License: BSD-3-Clause-LBNL
+
+
+# This is a script that analyses the simulation results from
+# the script `inputs_field_PEC_mr_3d`. This simulates a sinusoindal wave.
+# The electric field (Ey) is a standing wave due to the PEC boundary condition,
+# and as a result, the minimum and maximum value after reflection would be two times the value at initialization due to constructive interference.
+# Additionally, the value of Ey at the boundary must be equal to zero.
+import sys
+import re
+import matplotlib
+matplotlib.use('Agg')
+import matplotlib.pyplot as plt
+import yt
+yt.funcs.mylog.setLevel(50)
+import numpy as np
+from scipy.constants import e, m_e, epsilon_0, c
+sys.path.insert(1, '../../../../warpx/Regression/Checksum/')
+import checksumAPI
+
+# this will be the name of the plot file
+fn = sys.argv[1]
+
+# Parameters (these parameters must match the parameters in `inputs.multi.rt`)
+Ey_in = 1.e5
+E_th = 2.0*Ey_in
+
+# Read the file
+ds = yt.load(fn)
+
+t0 = ds.current_time.to_value()
+data = ds.covering_grid(level=0, left_edge=ds.domain_left_edge,
+                                    dims=ds.domain_dimensions)
+Ey_array = data['Ey'].to_ndarray()
+max_Ey_sim = Ey_array.max()
+min_Ey_sim = Ey_array.min()
+max_Ey_error = abs(max_Ey_sim-E_th)/abs(E_th)
+min_Ey_error = abs(min_Ey_sim - (-E_th))/abs(E_th)
+print('Ey_max is %s: Max Eth is %s : Max error for Ey_max: %.2e' %(max_Ey_sim,E_th,max_Ey_error))
+print('Ey_min is %s: Min Eth is %s : Max error for Ey_min: %.2e' %(min_Ey_sim,(-E_th), min_Ey_error))
+error_rel = 0.
+max_Ey_error_rel = max( error_rel, max_Ey_error )
+min_Ey_error_rel = max( error_rel, min_Ey_error )
+
+# Plot the last field from the loop (Ez at iteration 40)
+# Plot the last field from the loop (Ez at iteration 40)
+field = 'Ey'
+xCell = ds.domain_dimensions[0]
+yCell = ds.domain_dimensions[1]
+zCell = ds.domain_dimensions[2]
+z_array = data['z'].to_ndarray()
+
+plt.figure(figsize=(8,4))
+plt.plot(z_array[int(xCell/2),int(yCell/2),:]-z_array[int(xCell/2),int(yCell/2),int(zCell/2)],Ey_array[int(xCell/2),int(yCell/2),:])
+plt.ylim(-1.2e-3, 1.2e-3)
+plt.ylim(-2.2e5, 2.2e5)
+plt.xlim(-4.0e-6, 4.000001e-6)
+plt.xticks(np.arange(-4.e-6,4.000001e-6, step=1e-6))
+plt.xlabel('z (m)')
+plt.ylabel(field +' (V/m)')
+plt.tight_layout()
+plt.savefig('Ey_pec_analysis_mr.png')
+
+tolerance_rel = 0.05
+
+print("min_Ey_error_rel    : " + str(min_Ey_error_rel))
+print("max_Ey_error_rel    : " + str(max_Ey_error_rel))
+print("tolerance_rel: " + str(tolerance_rel))
+
+assert( max_Ey_error_rel < tolerance_rel )
+assert( min_Ey_error_rel < tolerance_rel )
+
+test_name = fn[:-9] # Could also be os.path.split(os.getcwd())[1]
+
+if re.search( 'single_precision', fn ):
+    checksumAPI.evaluate_checksum(test_name, fn, rtol=1.e-3)
+else:
+    checksumAPI.evaluate_checksum(test_name, fn)