diff options
author | 2021-11-19 00:31:18 -0800 | |
---|---|---|
committer | 2021-11-19 00:31:18 -0800 | |
commit | 45ef533c39a55d05afbb9872659963d19b79f463 (patch) | |
tree | 65d7756277fd82059e8267b92a65e16d7865590f /Examples/Modules | |
parent | ff8b73f58d71762e8d81be797c8afd1519d79c2a (diff) | |
download | WarpX-45ef533c39a55d05afbb9872659963d19b79f463.tar.gz WarpX-45ef533c39a55d05afbb9872659963d19b79f463.tar.zst WarpX-45ef533c39a55d05afbb9872659963d19b79f463.zip |
1D3V Cartesian Support (#2307)
* Build System: Add 1D Geometry
* test PR
* test PR
* 1D cartesian yee algorithm
* fixed typo
* Fixes for PML
* 1D support related multiple changes
* Fix compilation
* change 1D to 1D_Z
* 1D Field Gather and typo fix
* 1D Charge Deposition
* Particle Pusher
* multiple changes related to 1D
* 1D diagnostics and initialization
* PlasmaInjector and PEC fixes for 1D
* clean-up delete diags file
* mobility 1D laser particle and bilinear filter
* deleted diags files
* update laser particle weight formula
* delete diags files
* Azure: Add 1D Cartesian Runner
* 1D fixes for FieldProbe
* Update Docs/source/developers/dimensionality.rst
Co-authored-by: Remi Lehe <remi.lehe@normalesup.org>
* 1d laser injection and langmuir test input files
* 1d tests
* clean up : delete print statements
* analyse simulation result for laser injection and Langmuir tests
* EOL
* delete input files for which there are no automated tests
* delete input files for which there are no automated tests
* add ignore_unused to remove warnings
* remove space
* Fix compilation issues
* fix error : macro name must be an identifier
* Small bug fix
* cleanup Python script for analysis
* bug fix
* bug fix
* Update ParticleProbe: Check 1D in-domain
* Update Source/Make.WarpX
* Update .azure-pipelines.yml
* Add USE_OPENPMD=FALSE to .azure-pipeline.yml
* resolve conflict
* resolve conflict
* fix typo
* Correct out-of-bound access
* Fix Particle BC in WarpXParticleContainer and correct path to checksumAPI in python analysis scripts
* EOL
* Fix bug : accessing out of bound index of cell in 1D
* remove 1d test for cartesian3d
* Fix CI check
* Slight style change
* Address review comments
* Fix GPU compilation Filter.cpp
* Fix CI
* Fix Indentation
* Address review comments
* More consistent ifdef for dimension bigger than 1
* Update Examples/Tests/Langmuir/analysis_langmuir_multi_1d.py
Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
* Update GNUmakefile
Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
* Update Regression/prepare_file_ci.py
Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
* Update Source/FieldSolver/FiniteDifferenceSolver/FiniteDifferenceAlgorithms/CartesianCKCAlgorithm.H
Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
* Update Source/FieldSolver/FiniteDifferenceSolver/FiniteDifferenceAlgorithms/CartesianNodalAlgorithm.H
Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
* Update Source/Filter/Filter.cpp
Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
* Update Source/Filter/Filter.cpp
Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
* Update Source/Filter/Filter.cpp
Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
* Update Source/Filter/Filter.cpp
Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
* Update Source/Initialization/PlasmaInjector.cpp
Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
* Update Source/Initialization/PlasmaInjector.cpp
Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
* add comment inline to explain twice push_back
* Add amrex::Abort for NCIGodfreyFilter
Co-authored-by: Axel Huebl <axel.huebl@plasma.ninja>
Co-authored-by: Prabhat Kumar <prabhatkumar@kraken.dhcp.lbl.gov>
Co-authored-by: Remi Lehe <remi.lehe@normalesup.org>
Co-authored-by: Remi Lehe <rlehe@lbl.gov>
Diffstat (limited to 'Examples/Modules')
-rwxr-xr-x | Examples/Modules/laser_injection/analysis_1d.py | 175 | ||||
-rw-r--r-- | Examples/Modules/laser_injection/inputs_1d_rt | 62 |
2 files changed, 237 insertions, 0 deletions
diff --git a/Examples/Modules/laser_injection/analysis_1d.py b/Examples/Modules/laser_injection/analysis_1d.py new file mode 100755 index 000000000..6606938df --- /dev/null +++ b/Examples/Modules/laser_injection/analysis_1d.py @@ -0,0 +1,175 @@ +#! /usr/bin/env python + +# Copyright 2021 Prabhat Kumar, Remi Lehe +# +# This file is part of WarpX. +# +# License: BSD-3-Clause-LBNL + +# This file is part of the WarpX automated test suite. Its purpose is to test the +# injection of a Gaussian laser pulse from an antenna in a 1D simulation. +# The test calculates the envelope of each component of the laser pulse at the end of +# the simulation and it compares it with theory. It also checks that the +# central frequency of the Fourier transform is the expected one. + +import yt +import sys +import matplotlib +matplotlib.use('Agg') +import matplotlib.pyplot as plt +import numpy as np +from scipy.signal import hilbert +sys.path.insert(1, '../../../../warpx/Regression/Checksum/') +import checksumAPI + +# Maximum acceptable error for this test +relative_error_threshold = 0.05 + +# A small number +small_num = 1.0e-8 + +# Physical parameters +um = 1.e-6 +fs = 1.e-15 +c = 299792458 + +# Parameters of the gaussian beam +wavelength = 1.*um +w0 = 5.*um +tt = 10.*fs +t_c = 24.*fs +E_max = 4e12 + +# laser direction +dir_vector = np.array([0,0,1.0]) +dir_vector /= np.linalg.norm(dir_vector) + + +# polarization vector +pol_vector = np.array([1.0,1.0,0.0]) +pol_vector /= np.linalg.norm(pol_vector) + +# Calculates the envelope of a Gaussian beam +def gauss_env(T,Z): + '''Function to compute the theory for the envelope + ''' + inv_tau2 = 1./tt/tt + exp_arg = - inv_tau2 / c/c * (Z-T*c)*(Z-T*c) + return E_max * np.real(np.exp(exp_arg)) + +# Checks envelope and central frequency for a given laser component +def check_component(data, component, t_env_theory, coeff,Z,dz): + print("*** Checking " + component + " ***") + field = data['boxlib', component].v.squeeze() + env = abs(hilbert(field)) + + env_theory = t_env_theory*np.abs(coeff) + + # Plot results + fig = plt.figure(figsize=(12,6)) + + ax1 = fig.add_subplot(221) + ax1.set_title('PIC field') + ax1.plot(Z,field) + + ax2 = fig.add_subplot(222) + ax2.set_title('PIC envelope') + ax2.plot(Z,env) + + ax3 = fig.add_subplot(223) + ax3.set_title('Theory envelope') + ax3.plot(Z,env_theory, label="theory") + ax3.plot(Z,env, label="simulation") + ax3.legend(loc="upper right") + + ax4 = fig.add_subplot(224) + ax4.set_title('Difference') + ax4.plot(Z,env-env_theory) + + plt.tight_layout() + plt.savefig("plt_" + component + ".png", bbox_inches='tight') + + if(np.abs(coeff) < small_num): + is_field_zero = np.sum(np.abs(env)) < small_num + if is_field_zero : + print("[OK] Field component expected to be 0 is ~ 0") + else : + print("[FAIL] Field component expected to be 0 is NOT ~ 0") + assert(is_field_zero) + print("******\n") + return + + fft_field = np.fft.fft(field) + + freq_cols = np.fft.fftfreq(fft_field.shape[0],dz/c) + + pos_max = np.unravel_index(np.abs(fft_field).argmax(), fft_field.shape) + + freq = np.abs(freq_cols[pos_max[0]]) + exp_freq = c/wavelength + + relative_error_freq = np.abs(freq-exp_freq)/exp_freq + is_freq_ok = relative_error_freq < relative_error_threshold + if is_freq_ok : + print("[OK] Relative error frequency: {:6.3f} %".format(relative_error_freq*100)) + else : + print("[FAIL] Relative error frequency: {:6.3f} %".format(relative_error_freq*100)) + assert(is_freq_ok) + + print("******\n") + + relative_error_env = np.sum(np.abs(env-env_theory)) / np.sum(np.abs(env_theory)) + is_env_ok = relative_error_env < relative_error_threshold + if is_env_ok : + print("[OK] Relative error envelope: {:6.3f} %".format(relative_error_env*100)) + else : + print("[FAIL] Relative error envelope: {:6.3f} %".format(relative_error_env*100)) + assert(is_env_ok) + +def check_laser(filename): + ds = yt.load(filename) + + # yt 4.0+ has rounding issues with our domain data: + # RuntimeError: yt attempted to read outside the boundaries + # of a non-periodic domain along dimension 0. + if 'force_periodicity' in dir(ds): ds.force_periodicity() + + z = np.linspace( + ds.domain_left_edge[0].v, + ds.domain_right_edge[0].v, + ds.domain_dimensions[0]) + + dz = (ds.domain_right_edge[0].v-ds.domain_left_edge[0].v)/(ds.domain_dimensions[0]-1) + + # Compute the theory for envelope + env_theory = gauss_env(+t_c-ds.current_time.to_value(),z)+gauss_env(-t_c+ds.current_time.to_value(),z) + + # Read laser field in PIC simulation, and compute envelope + all_data_level_0 = ds.covering_grid(level=0, left_edge=ds.domain_left_edge, dims=ds.domain_dimensions) + + b_vector = np.cross(dir_vector, pol_vector) + + components = ["Ex", "Ey", "Ez", "Bx", "By", "Bz"] + coeffs = [ + pol_vector[0], + pol_vector[1], + pol_vector[2], + b_vector[0], + b_vector[1], + b_vector[2]] + + field_facts = [1, 1, 1, 1/c, 1/c, 1/c] + + for comp, coeff, field_fact in zip(components, coeffs, field_facts): + check_component(all_data_level_0, comp, field_fact*env_theory, coeff, z, dz) + +def main(): + filename_end = sys.argv[1] + + check_laser(filename_end) + + test_name = filename_end[:-9] # Could also be os.path.split(os.getcwd())[1] + checksumAPI.evaluate_checksum(test_name, filename_end) + +if __name__ == "__main__": + main() diff --git a/Examples/Modules/laser_injection/inputs_1d_rt b/Examples/Modules/laser_injection/inputs_1d_rt new file mode 100644 index 000000000..eb947039c --- /dev/null +++ b/Examples/Modules/laser_injection/inputs_1d_rt @@ -0,0 +1,62 @@ +# Maximum number of time steps +max_step = 240 + +# number of grid points +amr.n_cell = 352 + +# Maximum allowable size of each subdomain in the problem domain; +# this is used to decompose the domain for parallel calculations. +amr.max_grid_size = 32 + +# Maximum level in hierarchy (for now must be 0, i.e., one level in total) +amr.max_level = 0 + +# Geometry +geometry.coord_sys = 0 # 0: Cartesian +geometry.prob_lo = -15.e-6 # physical domain +geometry.prob_hi = 15.e-6 + +boundary.field_lo = pec +boundary.field_hi = pec + +warpx.serialize_ics = 1 + +# Verbosity +warpx.verbose = 1 + +# Algorithms +algo.current_deposition = esirkepov +warpx.use_filter = 0 + +# CFL +warpx.cfl = 0.9 + +# Order of particle shape factors +algo.particle_shape = 1 + +# Laser +lasers.names = laser1 +laser1.profile = Gaussian +laser1.position = 0.e-6 0.e-6 0.e-6 # This point is on the laser plane +laser1.direction = 0. 0. 1. # The plane normal direction +laser1.polarization = 1. 1. 0. # The main polarization vector +laser1.e_max = 4.e12 # Maximum amplitude of the laser field (in V/m) +laser1.wavelength = 1.0e-6 # The wavelength of the laser (in meters) +laser1.profile_waist = 5.e-6 # The waist of the laser (in meters) +laser1.profile_duration = 10.e-15 # The duration of the laser (in seconds) +laser1.profile_t_peak = 24.e-15 # The time at which the laser reaches its peak (in seconds) +laser1.profile_focal_distance = 13.109e-6 # Focal distance from the antenna (in meters) + # With this focal distance the laser is at focus + # at the end of the simulation. + +# Diagnostics +diagnostics.diags_names = diag1 +diag1.intervals = 20 +diag1.diag_type = Full + +# Moving window +warpx.do_moving_window = 1 +warpx.moving_window_dir = z +warpx.moving_window_v = 1.0 # in units of the speed of light +warpx.start_moving_window_step = 20 +warpx.end_moving_window_step = 200 |