ADI

.clangd

@@ -0,0 +1,46 @@
+# Match GNU Make flags for Bin/main3d.gnu.TPROF.MTMPI.ex (see `make -n …MacroscopicEvolveADI.o`).
+# Relative -I paths need compile_commands.json so clangd uses artemis/ as the
+# compile working directory (otherwise -I is resolved from each source file dir).
+CompileFlags:
+  CompilationDatabase: .
+  Add:
+    - -I.
+    - -I./Source
+    - -I./Source/Particles/Filter
+    - -Itmp_build_dir/s/3d.gnu.TPROF.MTMPI.EXE
+    - -I../amrex/Src/Base
+    - -I../amrex/Src/Base/Parser
+    - -I../amrex/Src/Particle
+    - -I../amrex/Src/Boundary
+    - -I../amrex/Src/AmrCore
+    - -I../amrex/Src/LinearSolvers/MLMG
+    - -I../amrex/Tools/C_scripts
+    - -I/usr/local/Cellar/open-mpi/5.0.8/include
+    - -DAMREX_EXPORT_DYNAMIC
+    - -DAMREX_TINY_PROFILING
+    - -DBL_USE_MPI
+    - -DAMREX_USE_MPI
+    - -DAMREX_MPI_THREAD_MULTIPLE
+    - -DBL_NO_FORT
+    - -DAMREX_GPU_MAX_THREADS=0
+    - -DBL_SPACEDIM=3
+    - -DAMREX_SPACEDIM=3
+    - -DBL_FORT_USE_UNDERSCORE
+    - -DAMREX_FORT_USE_UNDERSCORE
+    - -DBL_Darwin
+    - -DAMREX_Darwin
+    - -DAMREX_PARTICLES
+    - -DNDEBUG
+    - -DOMPI_SKIP_MPICXX
+    - -DWARPX_DIM_3D
+    - -DWARPX_ZINDEX=2
+    - -DPICSAR_NO_ASSUMED_ALIGNMENT
+    - -DWARPX
+
+---
+If:
+  PathMatch: .*\.(cpp|cc|cxx|hpp|ipp|hh|hxx|H)$
+CompileFlags:
+  Add:
+    - -xc++
+    - -std=c++20

.gitignore

@@ -12,6 +12,14 @@ AMReX_buildInfo.cpp
 Backtrace.*
 *.vth
 *.sw[opqrs]
+warpx_used_inputs
+diag*/
+wave_freespace_run/
+*_run/
+*.mp4
+*.mov
+*.avi
+*.gif
 
 ####################
 # AMReX data files #
@@ -70,3 +78,5 @@ cmake-build-*/
 .DS_Store
 .AppleDouble
 .LSOverride
+
+.cache

Examples/Tests/circuits/inputs_wave_freespace

@@ -1,10 +1,11 @@
 ################################
 ####### GENERAL PARAMETERS ######
 #################################
+# 3D mesh, 1D pulse along z (fine direction)
 max_step = 250
-amr.n_cell = 32 32 32
-amr.max_grid_size = 32
-amr.blocking_factor = 32
+amr.n_cell = 8 8 1024
+amr.max_grid_size = 1024
+amr.blocking_factor = 8
 geometry.dims = 3
 geometry.prob_lo     =  -4.e-6 -4.e-6 -4.e-6
 geometry.prob_hi     =   4.e-6  4.e-6  4.e-6
@@ -15,12 +16,14 @@ boundary.field_hi = periodic periodic periodic
 #################################
 ############ NUMERICS ###########
 #################################
-warpx.verbose = 1
+warpx.verbose = 0
 warpx.use_filter = 0
-warpx.cfl = 0.9
+warpx.cfl = 50
 
 algo.em_solver_medium = macroscopic # vacuum/macroscopic
 
+algo.time_stepping_scheme = adi # fdtd or adi
+
 algo.macroscopic_sigma_method = laxwendroff # laxwendroff or backwardeuler
 
 macroscopic.sigma_function(x,y,z) = "sigma_0"
@@ -42,15 +45,16 @@ my_constants.sigma_0 = 0.0
 my_constants.epsilon_0 = 8.8541878128e-12
 my_constants.mu_0 = 1.25663706212e-06
 
+# Ex / By packet propagating along +z
 warpx.E_ext_grid_init_style = parse_E_ext_grid_function
-warpx.Ex_external_grid_function(x,y,z) = "1.e5*exp(-(y)**2/(L/32)**2)*cos(2*pi*(y)/wavelength)"
+warpx.Ex_external_grid_function(x,y,z) = "1.e5*exp(-(z)**2/(L/32)**2)*cos(2*pi*(z)/wavelength)"
 warpx.Ey_external_grid_function(x,y,z) = "0."
 warpx.Ez_external_grid_function(x,y,z) = "0."
 
 warpx.B_ext_grid_init_style = parse_B_ext_grid_function
 warpx.Bx_external_grid_function(x,y,z) = "0."
-warpx.By_external_grid_function(x,y,z) = "0."
-warpx.Bz_external_grid_function(x,y,z) = "(-1.e5*exp(-(y)**2/(L/32)**2)*cos(2*pi*(y)/wavelength)/c)"
+warpx.By_external_grid_function(x,y,z) = "(1.e5*exp(-(z)**2/(L/32)**2)*cos(2*pi*(z)/wavelength)/c)"
+warpx.Bz_external_grid_function(x,y,z) = "0."
 
 warpx.E_excitation_on_grid_style = parse_E_excitation_grid_function
 

Examples/Tests/pec/inputs_field_PEC_3d

@@ -3,14 +3,14 @@
 # standing wave.
 
 # max step
-max_step = 125
+max_step = 1000
 
 # number of grid points
-amr.n_cell = 32 32 256
+amr.n_cell = 4 4 1024
 
 # Maximum allowable size of each subdomain
 amr.max_grid_size = 1024
-amr.blocking_factor = 32
+amr.blocking_factor = 4
 
 amr.max_level = 0
 
@@ -26,12 +26,21 @@ boundary.field_hi = periodic periodic pec
 warpx.serialize_initial_conditions = 1
 
 # Verbosity
-warpx.verbose = 1
+warpx.verbose = 0
 
 # Algorithms
 algo.current_deposition = esirkepov
+algo.em_solver_medium = macroscopic
+algo.time_stepping_scheme = adi
+algo.macroscopic_sigma_method = laxwendroff
+
+# Vacuum material properties for the macroscopic ADI field solve
+macroscopic.sigma_function(x,y,z) = "0.0"
+macroscopic.epsilon_function(x,y,z) = "8.8541878128e-12"
+macroscopic.mu_function(x,y,z) = "1.25663706212e-06"
+
 # CFL
-warpx.cfl = 0.9
+warpx.cfl = 2
 
 
 my_constants.z1 = -2.e-6
@@ -48,7 +57,7 @@ warpx.By_external_grid_function(x,y,z)= "0."
 warpx.Bz_external_grid_function(x,y,z) = "0."
 
 # Diagnostics
-diagnostics.diags_names = diag1
-diag1.intervals = 125
-diag1.diag_type = Full
-diag1.fields_to_plot = Ey Bx
+diagnostics.diags_names = plt
+plt.intervals = 10
+plt.diag_type = Full
+plt.fields_to_plot = Ex Ey Ez Bx By Bz

Source/Evolve/WarpXEvolve.cpp

@@ -150,7 +150,8 @@ WarpX::Evolve (int numsteps)
                 m_london->EvolveLondonJ(-0.5_rt*dt[0]); // J^(n) to J^(n-1/2) using E^(n)
                 FillBoundaryJ(guard_cells.ng_alloc_EB);
             }
-            if (use_lumped_inductor) {
+            if (use_lumped_inductor &&
+                macroscopic_time_integrator_algo != MacroscopicTimeSteppingScheme::ADI) {
                 m_inductor->EvolveInductorJ(-0.5_rt*dt[0]); // J^(n) to J^(n-1/2) using E^(n)
                 FillBoundaryJ(guard_cells.ng_alloc_EB);
             }
@@ -427,7 +428,8 @@ WarpX::OneStep_nosub (Real cur_time)
         // fill boundary here
     }
 #endif
-    if (use_lumped_inductor == 1) {
+    if (use_lumped_inductor == 1 &&
+        macroscopic_time_integrator_algo != MacroscopicTimeSteppingScheme::ADI) {
         m_inductor->EvolveInductorJ(dt[0]); // J^(n-1/2) to J^(n+1/2) using E^(n)
         FillBoundaryJ(guard_cells.ng_alloc_EB);
     }
@@ -486,11 +488,19 @@ WarpX::OneStep_nosub (Real cur_time)
         EvolveG(0.5_rt * dt[0], DtType::FirstHalf);
         FillBoundaryF(guard_cells.ng_FieldSolverF);
         FillBoundaryG(guard_cells.ng_FieldSolverG);
+        if (WarpX::em_solver_medium == MediumForEM::Macroscopic &&
+            WarpX::macroscopic_time_integrator_algo == MacroscopicTimeSteppingScheme::ADI) {
+#ifdef WARPX_MAG_LLG
+            WARPX_ABORT_WITH_MESSAGE("Macroscopic ADI is not implemented with WARPX_MAG_LLG.");
+#else
+            MacroscopicEvolveADI(dt[0]); // We now have E^{n+1} and B^{n+1}
+#endif
+        } else {
 #ifndef WARPX_MAG_LLG
         EvolveB(0.5_rt * dt[0], DtType::FirstHalf); // We now have B^{n+1/2}
         FillBoundaryB(guard_cells.ng_FieldSolver, WarpX::sync_nodal_points);
         // ApplyExternalFieldExcitation
-        ApplyExternalFieldExcitationOnGrid(ExternalFieldType::BfieldExternal, DtType::FirstHalf); // apply B external excitation; soft source to be fixed
+        ApplyExternalFieldExcitationOnGrid(ExternalFieldType::BfieldExternal, DtType::FirstHalf); // B soft/hard at t^{n+1/2}
 #endif
 
 #ifdef WARPX_MAG_LLG
@@ -550,8 +560,9 @@ WarpX::OneStep_nosub (Real cur_time)
 #ifndef WARPX_MAG_LLG
         EvolveB(0.5_rt * dt[0], DtType::SecondHalf); // We now have B^{n+1}
         // ApplyExternalFieldExcitation
-        ApplyExternalFieldExcitationOnGrid(ExternalFieldType::BfieldExternal, DtType::SecondHalf); // redundant for hs; need to fix the way to increment ss
+        ApplyExternalFieldExcitationOnGrid(ExternalFieldType::BfieldExternal, DtType::SecondHalf); // B soft/hard at t^{n+1}
 #endif
+        }
 
 #ifdef WARPX_MAG_LLG
 #ifndef WARPX_DIM_RZ

Source/FieldSolver/FiniteDifferenceSolver/CMakeLists.txt

@@ -10,6 +10,7 @@ target_sources(WarpX
     EvolveG.cpp
     EvolveECTRho.cpp
     FiniteDifferenceSolver.cpp
+    MacroscopicEvolveADI.cpp
     MacroscopicEvolveE.cpp
     MacroscopicEvolveEPML.cpp
     ApplySilverMuellerBoundary.cpp

Source/FieldSolver/FiniteDifferenceSolver/FiniteDifferenceSolver.H

@@ -18,6 +18,7 @@
 #include <AMReX_REAL.H>
 
 #include <AMReX_BaseFwd.H>
+#include <AMReX_Periodicity.H>
 
 #include <array>
 #include <memory>
@@ -127,6 +128,21 @@ class FiniteDifferenceSolver
                             amrex::Real const dt,
                             std::unique_ptr<MacroscopicProperties> const& macroscopic_properties);
 #ifndef WARPX_DIM_RZ
+#ifndef WARPX_MAG_LLG
+        /**
+          * \brief Periodic-only ADI update for macroscopic freespace cases.
+          */
+        void MacroscopicEvolveADI ( std::array< std::unique_ptr<amrex::MultiFab>, 3>& Efield,
+                            std::array< std::unique_ptr<amrex::MultiFab>, 3>& Bfield,
+                            std::array< std::array<std::unique_ptr<amrex::MultiFab>, 3>, 3>& Efield_adi,
+                            std::array< std::array<std::unique_ptr<amrex::MultiFab>, 3>, 3>& Bfield_adi,
+                            std::array< std::array<std::unique_ptr<amrex::MultiFab>, 3>, 3> const& PEC_adi,
+                            amrex::Real const dt,
+                            amrex::Periodicity const& periodicity,
+                            std::unique_ptr<MacroscopicProperties> const& macroscopic_properties);
+#endif
+#endif
+#ifndef WARPX_DIM_RZ
 #ifdef WARPX_MAG_LLG
         /**
           * \brief Macroscopic M-update, H-update and B=mu_o(H+M) computation for non-vacuum medium using finite-difference algorithm