diff --git a/.github/workflows/c-cpp.yml b/.github/workflows/c-cpp.yml index c7c66e71..b300b81a 100644 --- a/.github/workflows/c-cpp.yml +++ b/.github/workflows/c-cpp.yml @@ -28,6 +28,8 @@ jobs: - name: Run unit tests run: | cd build + ./unit_tests + OMP_NUM_THREADS=1 ./unit_tests -r junit > unit_tests_report.xml - name: Code coverage diff --git a/examples_hpc/lattice/lattice_S19_n10.cfg b/examples_hpc/lattice/lattice_S19_n10.cfg deleted file mode 100644 index 1cdc14f0..00000000 --- a/examples_hpc/lattice/lattice_S19_n10.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.1024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_19_n10 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_19_n10 -MESH_FILE = mesh/lattice_n10.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 19 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S19_n20.cfg b/examples_hpc/lattice/lattice_S19_n20.cfg deleted file mode 100644 index 8b115f7b..00000000 --- a/examples_hpc/lattice/lattice_S19_n20.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.2024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_19_n20 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_19_n20 -MESH_FILE = mesh/lattice_n20.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 19 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S19_n40.cfg b/examples_hpc/lattice/lattice_S19_n40.cfg deleted file mode 100644 index e168fd3a..00000000 --- a/examples_hpc/lattice/lattice_S19_n40.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.2024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_19_n40 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_19_n40 -MESH_FILE = mesh/lattice_n40.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 19 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S19_n80.cfg b/examples_hpc/lattice/lattice_S19_n80.cfg deleted file mode 100644 index d6afd57a..00000000 --- a/examples_hpc/lattice/lattice_S19_n80.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.8024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_19_n80 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_19_n80 -MESH_FILE = mesh/lattice_n80.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 19 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S29_n20.cfg b/examples_hpc/lattice/lattice_S29_n20.cfg deleted file mode 100644 index 57ba1dea..00000000 --- a/examples_hpc/lattice/lattice_S29_n20.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.2024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_29_n20 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_29_n20 -MESH_FILE = mesh/lattice_n20.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 29 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S29_n40.cfg b/examples_hpc/lattice/lattice_S29_n40.cfg deleted file mode 100644 index 5bc3ce55..00000000 --- a/examples_hpc/lattice/lattice_S29_n40.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.4024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_29_n40 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_29_n40 -MESH_FILE = mesh/lattice_n40.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 29 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S29_n80.cfg b/examples_hpc/lattice/lattice_S29_n80.cfg deleted file mode 100644 index 756fcfa8..00000000 --- a/examples_hpc/lattice/lattice_S29_n80.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.8024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_29_n80 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_29_n80 -MESH_FILE = mesh/lattice_n80.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 29 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S3_n20.cfg b/examples_hpc/lattice/lattice_S3_n20.cfg deleted file mode 100644 index 169d6f20..00000000 --- a/examples_hpc/lattice/lattice_S3_n20.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.2024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_3_n20 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_3_n20 -MESH_FILE = mesh/lattice_n20.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 3 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S41_n20.cfg b/examples_hpc/lattice/lattice_S41_n20.cfg deleted file mode 100644 index 2829042f..00000000 --- a/examples_hpc/lattice/lattice_S41_n20.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.2024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_41_n20 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_41_n20 -MESH_FILE = mesh/lattice_n20.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 41 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S41_n40.cfg b/examples_hpc/lattice/lattice_S41_n40.cfg deleted file mode 100644 index 258e5c5f..00000000 --- a/examples_hpc/lattice/lattice_S41_n40.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.4024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_41_n40 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_41_n40 -MESH_FILE = mesh/lattice_n40.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 41 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S41_n80.cfg b/examples_hpc/lattice/lattice_S41_n80.cfg deleted file mode 100644 index dbb614b6..00000000 --- a/examples_hpc/lattice/lattice_S41_n80.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.8024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_41_n80 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_41_n80 -MESH_FILE = mesh/lattice_n80.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 41 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S53_n20.cfg b/examples_hpc/lattice/lattice_S53_n20.cfg deleted file mode 100644 index 3b91b11e..00000000 --- a/examples_hpc/lattice/lattice_S53_n20.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.2024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_53_n20 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_53_n20 -MESH_FILE = mesh/lattice_n20.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 53 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S53_n40.cfg b/examples_hpc/lattice/lattice_S53_n40.cfg deleted file mode 100644 index d2d4ad5d..00000000 --- a/examples_hpc/lattice/lattice_S53_n40.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.4024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_53_n40 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_53_n40 -MESH_FILE = mesh/lattice_n40.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 53 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S53_n80.cfg b/examples_hpc/lattice/lattice_S53_n80.cfg deleted file mode 100644 index 4cf390dc..00000000 --- a/examples_hpc/lattice/lattice_S53_n80.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.8024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_53_n80 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_53_n80 -MESH_FILE = mesh/lattice_n80.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 53 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S65_n20.cfg b/examples_hpc/lattice/lattice_S65_n20.cfg deleted file mode 100644 index b3613d0e..00000000 --- a/examples_hpc/lattice/lattice_S65_n20.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.2024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_65_n20 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_65_n20 -MESH_FILE = mesh/lattice_n20.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 65 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S65_n40.cfg b/examples_hpc/lattice/lattice_S65_n40.cfg deleted file mode 100644 index aa2645b4..00000000 --- a/examples_hpc/lattice/lattice_S65_n40.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.4024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_65_n40 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_65_n40 -MESH_FILE = mesh/lattice_n40.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 65 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S65_n80.cfg b/examples_hpc/lattice/lattice_S65_n80.cfg deleted file mode 100644 index adc3cf68..00000000 --- a/examples_hpc/lattice/lattice_S65_n80.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.8024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_65_n80 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_65_n80 -MESH_FILE = mesh/lattice_n80.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 65 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/lattice_S9_n20.cfg b/examples_hpc/lattice/lattice_S9_n20.cfg deleted file mode 100644 index 895b733f..00000000 --- a/examples_hpc/lattice/lattice_S9_n20.cfg +++ /dev/null @@ -1,48 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% Lattice Benchmarking File SN % -% Author % -% Date 08.01.2024 % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% -% ----IO specification ---- -% -OUTPUT_DIR = result -OUTPUT_FILE = lattice_quad_order_9_n20 -LOG_DIR = result/logs -LOG_FILE = lattice_quad_order_9_n20 -MESH_FILE = mesh/lattice_n20.su2 -% -% --- Problem definition --- -% -PROBLEM = LATTICE -TIME_FINAL = 2.8 -SPATIAL_DIM = 3 -SOURCE_MAGNITUDE = 1.0 -% -% ---- Solver specifications ---- -% -% Solver type -SOLVER = SN_SOLVER -% CFL number -CFL_NUMBER = 0.9 -% Reconstruction order -RECONS_ORDER = 2 -% -% ---- Boundary Conditions ---- -% -BC_NEUMANN = ( void ) -% -% ----- Quadrature Specification --- -% -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 9 -% -% ----- Output ---- -% -VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 0 -SCREEN_OUTPUT = (ITER, MASS,RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION ) -SCREEN_OUTPUT_FREQUENCY = 100 -HISTORY_OUTPUT = (ITER, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, CUR_PARTICLE_ABSORPTION, TOTAL_PARTICLE_ABSORPTION, MAX_PARTICLE_ABSORPTION) -HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/examples_hpc/lattice/mesh/lattice_rectangular_n20.geo b/examples_hpc/lattice/mesh/lattice_rectangular_n20.geo deleted file mode 100644 index 9721aecd..00000000 --- a/examples_hpc/lattice/mesh/lattice_rectangular_n20.geo +++ /dev/null @@ -1,1193 +0,0 @@ -cl_fine = 0.1; -n_recombine = 20; -n_prog = 1.1; -Point(1) = {-3.5, -3.5, 0, cl_fine}; -Point(2) = {3.5, -3.5, 0, cl_fine}; -Point(3) = {-3.5, 3.5, 0, cl_fine}; -Point(4) = {3.5, 3.5, 0, cl_fine}; - -// Inner grid -Point(5) = { -2.5, -2.5, 0, cl_fine }; -Point(6) = { -2.5, -2.0, 0, cl_fine }; -Point(7) = { -2.5, -1.5, 0, cl_fine }; -Point(8) = { -2.5, -1.0, 0, cl_fine }; -Point(9) = { -2.5, -0.5, 0, cl_fine }; -Point(10) = { -2.5, 0.0, 0, cl_fine }; -Point(11) = { -2.5, 0.5, 0, cl_fine }; -Point(12) = { -2.5, 1.0, 0, cl_fine }; -Point(13) = { -2.5, 1.5, 0, cl_fine }; -Point(14) = { -2.5, 2.0, 0, cl_fine }; -Point(15) = { -2.5, 2.5, 0, cl_fine }; -Point(16) = { -2.0, -2.5, 0, cl_fine }; -Point(17) = { -2.0, -2.0, 0, cl_fine }; -Point(18) = { -2.0, -1.5, 0, cl_fine }; -Point(19) = { -2.0, -1.0, 0, cl_fine }; -Point(20) = { -2.0, -0.5, 0, cl_fine }; -Point(21) = { -2.0, 0.0, 0, cl_fine }; -Point(22) = { -2.0, 0.5, 0, cl_fine }; -Point(23) = { -2.0, 1.0, 0, cl_fine }; -Point(24) = { -2.0, 1.5, 0, cl_fine }; -Point(25) = { -2.0, 2.0, 0, cl_fine }; -Point(26) = { -2.0, 2.5, 0, cl_fine }; -Point(27) = { -1.5, -2.5, 0, cl_fine }; -Point(28) = { -1.5, -2.0, 0, cl_fine }; -Point(29) = { -1.5, -1.5, 0, cl_fine }; -Point(30) = { -1.5, -1.0, 0, cl_fine }; -Point(31) = { -1.5, -0.5, 0, cl_fine }; -Point(32) = { -1.5, 0.0, 0, cl_fine }; -Point(33) = { -1.5, 0.5, 0, cl_fine }; -Point(34) = { -1.5, 1.0, 0, cl_fine }; -Point(35) = { -1.5, 1.5, 0, cl_fine }; -Point(36) = { -1.5, 2.0, 0, cl_fine }; -Point(37) = { -1.5, 2.5, 0, cl_fine }; -Point(38) = { -1.0, -2.5, 0, cl_fine }; -Point(39) = { -1.0, -2.0, 0, cl_fine }; -Point(40) = { -1.0, -1.5, 0, cl_fine }; -Point(41) = { -1.0, -1.0, 0, cl_fine }; -Point(42) = { -1.0, -0.5, 0, cl_fine }; -Point(43) = { -1.0, 0.0, 0, cl_fine }; -Point(44) = { -1.0, 0.5, 0, cl_fine }; -Point(45) = { -1.0, 1.0, 0, cl_fine }; -Point(46) = { -1.0, 1.5, 0, cl_fine }; -Point(47) = { -1.0, 2.0, 0, cl_fine }; -Point(48) = { -1.0, 2.5, 0, cl_fine }; -Point(49) = { -0.5, -2.5, 0, cl_fine }; -Point(50) = { -0.5, -2.0, 0, cl_fine }; -Point(51) = { -0.5, -1.5, 0, cl_fine }; -Point(52) = { -0.5, -1.0, 0, cl_fine }; -Point(53) = { -0.5, -0.5, 0, cl_fine }; -Point(54) = { -0.5, 0.0, 0, cl_fine }; -Point(55) = { -0.5, 0.5, 0, cl_fine }; -Point(56) = { -0.5, 1.0, 0, cl_fine }; -Point(57) = { -0.5, 1.5, 0, cl_fine }; -Point(58) = { -0.5, 2.0, 0, cl_fine }; -Point(59) = { -0.5, 2.5, 0, cl_fine }; -Point(60) = { 0.0, -2.5, 0, cl_fine }; -Point(61) = { 0.0, -2.0, 0, cl_fine }; -Point(62) = { 0.0, -1.5, 0, cl_fine }; -Point(63) = { 0.0, -1.0, 0, cl_fine }; -Point(64) = { 0.0, -0.5, 0, cl_fine }; -Point(65) = { 0.0, 0.0, 0, cl_fine }; -Point(66) = { 0.0, 0.5, 0, cl_fine }; -Point(67) = { 0.0, 1.0, 0, cl_fine }; -Point(68) = { 0.0, 1.5, 0, cl_fine }; -Point(69) = { 0.0, 2.0, 0, cl_fine }; -Point(70) = { 0.0, 2.5, 0, cl_fine }; -Point(71) = { 0.5, -2.5, 0, cl_fine }; -Point(72) = { 0.5, -2.0, 0, cl_fine }; -Point(73) = { 0.5, -1.5, 0, cl_fine }; -Point(74) = { 0.5, -1.0, 0, cl_fine }; -Point(75) = { 0.5, -0.5, 0, cl_fine }; -Point(76) = { 0.5, 0.0, 0, cl_fine }; -Point(77) = { 0.5, 0.5, 0, cl_fine }; -Point(78) = { 0.5, 1.0, 0, cl_fine }; -Point(79) = { 0.5, 1.5, 0, cl_fine }; -Point(80) = { 0.5, 2.0, 0, cl_fine }; -Point(81) = { 0.5, 2.5, 0, cl_fine }; -Point(82) = { 1.0, -2.5, 0, cl_fine }; -Point(83) = { 1.0, -2.0, 0, cl_fine }; -Point(84) = { 1.0, -1.5, 0, cl_fine }; -Point(85) = { 1.0, -1.0, 0, cl_fine }; -Point(86) = { 1.0, -0.5, 0, cl_fine }; -Point(87) = { 1.0, 0.0, 0, cl_fine }; -Point(88) = { 1.0, 0.5, 0, cl_fine }; -Point(89) = { 1.0, 1.0, 0, cl_fine }; -Point(90) = { 1.0, 1.5, 0, cl_fine }; -Point(91) = { 1.0, 2.0, 0, cl_fine }; -Point(92) = { 1.0, 2.5, 0, cl_fine }; -Point(93) = { 1.5, -2.5, 0, cl_fine }; -Point(94) = { 1.5, -2.0, 0, cl_fine }; -Point(95) = { 1.5, -1.5, 0, cl_fine }; -Point(96) = { 1.5, -1.0, 0, cl_fine }; -Point(97) = { 1.5, -0.5, 0, cl_fine }; -Point(98) = { 1.5, 0.0, 0, cl_fine }; -Point(99) = { 1.5, 0.5, 0, cl_fine }; -Point(100) = { 1.5, 1.0, 0, cl_fine }; -Point(101) = { 1.5, 1.5, 0, cl_fine }; -Point(102) = { 1.5, 2.0, 0, cl_fine }; -Point(103) = { 1.5, 2.5, 0, cl_fine }; -Point(104) = { 2.0, -2.5, 0, cl_fine }; -Point(105) = { 2.0, -2.0, 0, cl_fine }; -Point(106) = { 2.0, -1.5, 0, cl_fine }; -Point(107) = { 2.0, -1.0, 0, cl_fine }; -Point(108) = { 2.0, -0.5, 0, cl_fine }; -Point(109) = { 2.0, 0.0, 0, cl_fine }; -Point(110) = { 2.0, 0.5, 0, cl_fine }; -Point(111) = { 2.0, 1.0, 0, cl_fine }; -Point(112) = { 2.0, 1.5, 0, cl_fine }; -Point(113) = { 2.0, 2.0, 0, cl_fine }; -Point(114) = { 2.0, 2.5, 0, cl_fine }; -Point(115) = { 2.5, -2.5, 0, cl_fine }; -Point(116) = { 2.5, -2.0, 0, cl_fine }; -Point(117) = { 2.5, -1.5, 0, cl_fine }; -Point(118) = { 2.5, -1.0, 0, cl_fine }; -Point(119) = { 2.5, -0.5, 0, cl_fine }; -Point(120) = { 2.5, 0.0, 0, cl_fine }; -Point(121) = { 2.5, 0.5, 0, cl_fine }; -Point(122) = { 2.5, 1.0, 0, cl_fine }; -Point(123) = { 2.5, 1.5, 0, cl_fine }; -Point(124) = { 2.5, 2.0, 0, cl_fine }; -Point(125) = { 2.5, 2.5, 0, cl_fine }; - -// helper boundary points -Point(126) = { -2.5, -3.5, 0, cl_fine }; -Point(127) = { -2.0, -3.5, 0, cl_fine }; -Point(128) = { -1.5, -3.5, 0, cl_fine }; -Point(129) = { -1.0, -3.5, 0, cl_fine }; -Point(130) = { -0.5, -3.5, 0, cl_fine }; -Point(131) = { 0.0, -3.5, 0, cl_fine }; -Point(132) = { 0.5, -3.5, 0, cl_fine }; -Point(133) = { 1.0, -3.5, 0, cl_fine }; -Point(134) = { 1.5, -3.5, 0, cl_fine }; -Point(135) = { 2.0, -3.5, 0, cl_fine }; -Point(136) = { 2.5, -3.5, 0, cl_fine }; -Point(137) = { -2.5, 3.5, 0, cl_fine }; -Point(138) = { -2.0, 3.5, 0, cl_fine }; -Point(139) = { -1.5, 3.5, 0, cl_fine }; -Point(140) = { -1.0, 3.5, 0, cl_fine }; -Point(141) = { -0.5, 3.5, 0, cl_fine }; -Point(142) = { 0.0, 3.5, 0, cl_fine }; -Point(143) = { 0.5, 3.5, 0, cl_fine }; -Point(144) = { 1.0, 3.5, 0, cl_fine }; -Point(145) = { 1.5, 3.5, 0, cl_fine }; -Point(146) = { 2.0, 3.5, 0, cl_fine }; -Point(147) = { 2.5, 3.5, 0, cl_fine }; -Point(148) = { -3.5, -2.5, 0, cl_fine }; -Point(149) = { -3.5, -2.0, 0, cl_fine }; -Point(150) = { -3.5, -1.5, 0, cl_fine }; -Point(151) = { -3.5, -1.0, 0, cl_fine }; -Point(152) = { -3.5, -0.5, 0, cl_fine }; -Point(153) = { -3.5, 0.0, 0, cl_fine }; -Point(154) = { -3.5, 0.5, 0, cl_fine }; -Point(155) = { -3.5, 1.0, 0, cl_fine }; -Point(156) = { -3.5, 1.5, 0, cl_fine }; -Point(157) = { -3.5, 2.0, 0, cl_fine }; -Point(158) = { -3.5, 2.5, 0, cl_fine }; -Point(159) = { 3.5, -2.5, 0, cl_fine }; -Point(160) = { 3.5, -2.0, 0, cl_fine }; -Point(161) = { 3.5, -1.5, 0, cl_fine }; -Point(162) = { 3.5, -1.0, 0, cl_fine }; -Point(163) = { 3.5, -0.5, 0, cl_fine }; -Point(164) = { 3.5, 0.0, 0, cl_fine }; -Point(165) = { 3.5, 0.5, 0, cl_fine }; -Point(166) = { 3.5, 1.0, 0, cl_fine }; -Point(167) = { 3.5, 1.5, 0, cl_fine }; -Point(168) = { 3.5, 2.0, 0, cl_fine }; -Point(169) = { 3.5, 2.5, 0, cl_fine }; - -// Horizontal and vertical lines in inner grid -Line(5) = { 5, 16}; -Line(6) = { 6, 17}; -Line(7) = { 7, 18}; -Line(8) = { 8, 19}; -Line(9) = { 9, 20}; -Line(10) = { 10, 21}; -Line(11) = { 11, 22}; -Line(12) = { 12, 23}; -Line(13) = { 13, 24}; -Line(14) = { 14, 25}; -Line(15) = { 15, 26}; -Line(16) = { 27, 16}; -Line(17) = { 28, 17}; -Line(18) = { 29, 18}; -Line(19) = { 30, 19}; -Line(20) = { 31, 20}; -Line(21) = { 32, 21}; -Line(22) = { 33, 22}; -Line(23) = { 34, 23}; -Line(24) = { 35, 24}; -Line(25) = { 36, 25}; -Line(26) = { 37, 26}; -Line(27) = { 27, 38}; -Line(28) = { 28, 39}; -Line(29) = { 29, 40}; -Line(30) = { 30, 41}; -Line(31) = { 31, 42}; -Line(32) = { 32, 43}; -Line(33) = { 33, 44}; -Line(34) = { 34, 45}; -Line(35) = { 35, 46}; -Line(36) = { 36, 47}; -Line(37) = { 37, 48}; -Line(38) = { 49, 38}; -Line(39) = { 50, 39}; -Line(40) = { 51, 40}; -Line(41) = { 52, 41}; -Line(42) = { 53, 42}; -Line(43) = { 54, 43}; -Line(44) = { 55, 44}; -Line(45) = { 56, 45}; -Line(46) = { 57, 46}; -Line(47) = { 58, 47}; -Line(48) = { 59, 48}; -Line(49) = { 49, 60}; -Line(50) = { 50, 61}; -Line(51) = { 51, 62}; -Line(52) = { 52, 63}; -Line(53) = { 53, 64}; -Line(54) = { 54, 65}; -Line(55) = { 55, 66}; -Line(56) = { 56, 67}; -Line(57) = { 57, 68}; -Line(58) = { 58, 69}; -Line(59) = { 59, 70}; -Line(60) = { 71, 60}; -Line(61) = { 72, 61}; -Line(62) = { 73, 62}; -Line(63) = { 74, 63}; -Line(64) = { 75, 64}; -Line(65) = { 76, 65}; -Line(66) = { 77, 66}; -Line(67) = { 78, 67}; -Line(68) = { 79, 68}; -Line(69) = { 80, 69}; -Line(70) = { 81, 70}; -Line(71) = { 71, 82}; -Line(72) = { 72, 83}; -Line(73) = { 73, 84}; -Line(74) = { 74, 85}; -Line(75) = { 75, 86}; -Line(76) = { 76, 87}; -Line(77) = { 77, 88}; -Line(78) = { 78, 89}; -Line(79) = { 79, 90}; -Line(80) = { 80, 91}; -Line(81) = { 81, 92}; -Line(82) = { 93, 82}; -Line(83) = { 94, 83}; -Line(84) = { 95, 84}; -Line(85) = { 96, 85}; -Line(86) = { 97, 86}; -Line(87) = { 98, 87}; -Line(88) = { 99, 88}; -Line(89) = { 100, 89}; -Line(90) = { 101, 90}; -Line(91) = { 102, 91}; -Line(92) = { 103, 92}; -Line(93) = { 93, 104}; -Line(94) = { 94, 105}; -Line(95) = { 95, 106}; -Line(96) = { 96, 107}; -Line(97) = { 97, 108}; -Line(98) = { 98, 109}; -Line(99) = { 99, 110}; -Line(100) = { 100, 111}; -Line(101) = { 101, 112}; -Line(102) = { 102, 113}; -Line(103) = { 103, 114}; -Line(104) = { 115, 104}; -Line(105) = { 116, 105}; -Line(106) = { 117, 106}; -Line(107) = { 118, 107}; -Line(108) = { 119, 108}; -Line(109) = { 120, 109}; -Line(110) = { 121, 110}; -Line(111) = { 122, 111}; -Line(112) = { 123, 112}; -Line(113) = { 124, 113}; -Line(114) = { 125, 114}; -Line(126) = { 5, 6}; -Line(127) = { 7, 6}; -Line(128) = { 7, 8}; -Line(129) = { 9, 8}; -Line(130) = { 9, 10}; -Line(131) = { 11, 10}; -Line(132) = { 11, 12}; -Line(133) = { 13, 12}; -Line(134) = { 13, 14}; -Line(135) = { 15, 14}; -Line(136) = { 16, 17}; -Line(137) = { 18, 17}; -Line(138) = { 18, 19}; -Line(139) = { 20, 19}; -Line(140) = { 20, 21}; -Line(141) = { 22, 21}; -Line(142) = { 22, 23}; -Line(143) = { 24, 23}; -Line(144) = { 24, 25}; -Line(145) = { 26, 25}; - -Line(146) = { 27, 28}; -Line(147) = { 29, 28}; -Line(148) = { 29, 30}; -Line(149) = { 31, 30}; -Line(150) = { 31, 32}; -Line(151) = { 33, 32}; -Line(152) = { 33, 34}; -Line(153) = { 35, 34}; -Line(154) = { 35, 36}; -Line(155) = { 37, 36}; - -Line(156) = { 38, 39}; -Line(157) = { 40, 39}; -Line(158) = { 40, 41}; -Line(159) = { 42, 41}; -Line(160) = { 42, 43}; -Line(161) = { 44, 43}; -Line(162) = { 44, 45}; -Line(163) = { 46, 45}; -Line(164) = { 46, 47}; -Line(165) = { 48, 47}; - -Line(166) = { 49, 50}; -Line(167) = { 51, 50}; -Line(168) = { 51, 52}; -Line(169) = { 53, 52}; -Line(170) = { 53, 54}; -Line(171) = { 55, 54}; -Line(172) = { 55, 56}; -Line(173) = { 57, 56}; -Line(174) = { 57, 58}; -Line(175) = { 59, 58}; - -Line(176) = { 60, 61}; -Line(177) = { 62, 61}; -Line(178) = { 62, 63}; -Line(179) = { 64, 63}; -Line(180) = { 64, 65}; -Line(181) = { 66, 65}; -Line(182) = { 66, 67}; -Line(183) = { 68, 67}; -Line(184) = { 68, 69}; -Line(185) = { 70, 69}; - -Line(186) = { 71, 72}; -Line(187) = { 73, 72}; -Line(188) = { 73, 74}; -Line(189) = { 75, 74}; -Line(190) = { 75, 76}; -Line(191) = { 77, 76}; -Line(192) = { 77, 78}; -Line(193) = { 79, 78}; -Line(194) = { 79, 80}; -Line(195) = { 81, 80}; - -Line(196) = { 82, 83}; -Line(197) = { 84, 83}; -Line(198) = { 84, 85}; -Line(199) = { 86, 85}; -Line(200) = { 86, 87}; -Line(201) = { 88, 87}; -Line(202) = { 88, 89}; -Line(203) = { 90, 89}; -Line(204) = { 90, 91}; -Line(205) = { 92, 91}; - -Line(206) = { 93, 94}; -Line(207) = { 95, 94}; -Line(208) = { 95, 96}; -Line(209) = { 97, 96}; -Line(210) = { 97, 98}; -Line(211) = { 99, 98}; -Line(212) = { 99, 100}; -Line(213) = { 101, 100}; -Line(214) = { 101, 102}; -Line(215) = { 103, 102}; - - -Line(216) = { 104, 105}; -Line(217) = { 106, 105}; -Line(218) = { 106, 107}; -Line(219) = { 108, 107}; -Line(220) = { 108, 109}; -Line(221) = { 110, 109}; -Line(222) = { 110, 111}; -Line(223) = { 112, 111}; -Line(224) = { 112, 113}; -Line(225) = { 114, 113}; - - -Line(226) = { 115, 116}; -Line(227) = { 117, 116}; -Line(228) = { 117, 118}; -Line(229) = { 119, 118}; -Line(230) = { 119, 120}; -Line(231) = { 121, 120}; -Line(232) = { 121, 122}; -Line(233) = { 123, 122}; -Line(234) = { 123, 124}; -Line(235) = { 125, 124}; - - -// Define inner surfaces -//+ -Curve Loop(1) = { 5, 136, -6, -126}; -Plane Surface(1) ={1}; -Curve Loop(2) = { 6, -137, -7, 127}; -Plane Surface(2) ={2}; -Curve Loop(3) = { 7, 138, -8, -128}; -Plane Surface(3) ={3}; -Curve Loop(4) = { 8, -139, -9, 129}; -Plane Surface(4) ={4}; -Curve Loop(5) = { 9, 140, -10, -130}; -Plane Surface(5) ={5}; -Curve Loop(6) = { 10, -141, -11, 131}; -Plane Surface(6) ={6}; -Curve Loop(7) = { 11, 142, -12, -132}; -Plane Surface(7) ={7}; -Curve Loop(8) = { 12, -143, -13, 133}; -Plane Surface(8) ={8}; -Curve Loop(9) = { 13, 144, -14, -134}; -Plane Surface(9) ={9}; -Curve Loop(10) = { 14, -145, -15, 135}; -Plane Surface(10) ={10}; -Curve Loop(11) = { -16, 146, 17, -136}; -Plane Surface(11) ={11}; -Curve Loop(12) = { -17, -147, 18, 137}; -Plane Surface(12) ={12}; -Curve Loop(13) = { -18, 148, 19, -138}; -Plane Surface(13) ={13}; -Curve Loop(14) = { -19, -149, 20, 139}; -Plane Surface(14) ={14}; -Curve Loop(15) = { -20, 150, 21, -140}; -Plane Surface(15) ={15}; -Curve Loop(16) = { -21, -151, 22, 141}; -Plane Surface(16) ={16}; -Curve Loop(17) = { -22, 152, 23, -142}; -Plane Surface(17) ={17}; -Curve Loop(18) = { -23, -153, 24, 143}; -Plane Surface(18) ={18}; -Curve Loop(19) = {-24, 154, 25, -144}; -Plane Surface(19) ={19}; -Curve Loop(20) = { -25, -155, 26, 145}; -Plane Surface(20) ={20}; -Curve Loop(21) = { 27, 156, -28, -146}; -Plane Surface(21) ={21}; -Curve Loop(22) = { 28, -157, -29, 147}; -Plane Surface(22) ={22}; -Curve Loop(23) = { 29, 158, -30, -148}; -Plane Surface(23) ={23}; -Curve Loop(24) = { 30, -159, -31, 149}; -Plane Surface(24) ={24}; -Curve Loop(25) = { 31, 160, -32, -150}; -Plane Surface(25) ={25}; -Curve Loop(26) = { 32, -161, -33, 151}; -Plane Surface(26) ={26}; -Curve Loop(27) = { 33, 162, -34, -152}; -Plane Surface(27) ={27}; -Curve Loop(28) = { 34, -163, -35, 153}; -Plane Surface(28) ={28}; -Curve Loop(29) = { 35, 164, -36, -154}; -Plane Surface(29) ={29}; -Curve Loop(30) = { 36, -165, -37, 155}; -Plane Surface(30) ={30}; -Curve Loop(31) = { -38, 166, 39, -156}; -Plane Surface(31) ={31}; -Curve Loop(32) = { -39, -167, 40, 157}; -Plane Surface(32) ={32}; -Curve Loop(33) = { -40, 168, 41, -158}; -Plane Surface(33) ={33}; -Curve Loop(34) = { -41, -169, 42, 159}; -Plane Surface(34) ={34}; -Curve Loop(35) = { -42, 170, 43, -160}; -Plane Surface(35) ={35}; -Curve Loop(36) = { -43, -171, 44, 161}; -Plane Surface(36) ={36}; -Curve Loop(37) = { -44, 172, 45, -162}; -Plane Surface(37) ={37}; -Curve Loop(38) = { -45, -173, 46, 163}; -Plane Surface(38) ={38}; -Curve Loop(39) = { -46, 174, 47, -164}; -Plane Surface(39) ={39}; -Curve Loop(40) = { -47, -175, 48, 165}; -Plane Surface(40) ={40}; -Curve Loop(41) = { 49, 176, -50, -166}; -Plane Surface(41) ={41}; -Curve Loop(42) = { 50, -177, -51, 167}; -Plane Surface(42) ={42}; -Curve Loop(43) = { 51, 178, -52, -168}; -Plane Surface(43) ={43}; -Curve Loop(44) = { 52, -179, -53, 169}; -Plane Surface(44) ={44}; -Curve Loop(45) = { 53, 180, -54, -170}; -Plane Surface(45) ={45}; -Curve Loop(46) = { 54, -181, -55, 171}; -Plane Surface(46) ={46}; -Curve Loop(47) = { 55, 182, -56, -172}; -Plane Surface(47) ={47}; -Curve Loop(48) = { 56, -183, -57, 173}; -Plane Surface(48) ={48}; -Curve Loop(49) = { 57, 184, -58, -174}; -Plane Surface(49) ={49}; -Curve Loop(50) = { 58, -185, -59, 175}; -Plane Surface(50) ={50}; -Curve Loop(51) = { -60, 186, 61, -176}; -Plane Surface(51) ={51}; -Curve Loop(52) = { -61, -187, 62, 177}; -Plane Surface(52) ={52}; -Curve Loop(53) = { -62, 188, 63, -178}; -Plane Surface(53) ={53}; -Curve Loop(54) = { -63, -189, 64, 179}; -Plane Surface(54) ={54}; -Curve Loop(55) = { -64, 190, 65, -180}; -Plane Surface(55) ={55}; -Curve Loop(56) = { -65, -191, 66, 181}; -Plane Surface(56) ={56}; -Curve Loop(57) = { -66, 192, 67, -182}; -Plane Surface(57) ={57}; -Curve Loop(58) = { -67, -193, 68, 183}; -Plane Surface(58) ={58}; -Curve Loop(59) = { -68, 194, 69, -184}; -Plane Surface(59) ={59}; -Curve Loop(60) = { -69, -195, 70, 185}; -Plane Surface(60) ={60}; -Curve Loop(61) = { 71, 196, -72, -186}; -Plane Surface(61) ={61}; -Curve Loop(62) = { 72, -197, -73, 187}; -Plane Surface(62) ={62}; -Curve Loop(63) = { 73, 198, -74, -188}; -Plane Surface(63) ={63}; -Curve Loop(64) = { 74, -199, -75, 189}; -Plane Surface(64) ={64}; -Curve Loop(65) = { 75, 200, -76, -190}; -Plane Surface(65) ={65}; -Curve Loop(66) = { 76, -201, -77, 191}; -Plane Surface(66) ={66}; -Curve Loop(67) = { 77, 202, -78, -192}; -Plane Surface(67) ={67}; -Curve Loop(68) = { 78, -203, -79, 193}; -Plane Surface(68) ={68}; -Curve Loop(69) = { 79, 204, -80, -194}; -Plane Surface(69) ={69}; -Curve Loop(70) = { 80, -205, -81, 195}; -Plane Surface(70) ={70}; -Curve Loop(71) = { -82, 206, 83, -196}; -Plane Surface(71) ={71}; -Curve Loop(72) = { -83, -207, 84, 197}; -Plane Surface(72) ={72}; -Curve Loop(73) = { -84, 208, 85, -198}; -Plane Surface(73) ={73}; -Curve Loop(74) = { -85, -209, 86, 199}; -Plane Surface(74) ={74}; -Curve Loop(75) = { -86, 210, 87, -200}; -Plane Surface(75) ={75}; -Curve Loop(76) = { -87, -211, 88, 201}; -Plane Surface(76) ={76}; -Curve Loop(77) = { -88, 212, 89, -202}; -Plane Surface(77) ={77}; -Curve Loop(78) = { -89, -213, 90, 203}; -Plane Surface(78) ={78}; -Curve Loop(79) = { -90, 214, 91, -204}; -Plane Surface(79) ={79}; -Curve Loop(80) = { -91, -215, 92, 205}; -Plane Surface(80) ={80}; -Curve Loop(81) = { 93, 216, -94, -206}; -Plane Surface(81) ={81}; -Curve Loop(82) = { 94, -217, -95, 207}; -Plane Surface(82) ={82}; -Curve Loop(83) = { 95, 218, -96, -208}; -Plane Surface(83) ={83}; -Curve Loop(84) = { 96, -219, -97, 209}; -Plane Surface(84) ={84}; -Curve Loop(85) = { 97, 220, -98, -210}; -Plane Surface(85) ={85}; -Curve Loop(86) = { 98, -221, -99, 211}; -Plane Surface(86) ={86}; -Curve Loop(87) = { 99, 222, -100, -212}; -Plane Surface(87) ={87}; -Curve Loop(88) = { 100, -223, -101, 213}; -Plane Surface(88) ={88}; -Curve Loop(89) = { 101, 224, -102, -214}; -Plane Surface(89) ={89}; -Curve Loop(90) = { 102, -225, -103, 215}; -Plane Surface(90) ={90}; -Curve Loop(91) = { -104, 226, 105, -216}; -Plane Surface(91) ={91}; -Curve Loop(92) = { -105, -227, 106, 217}; -Plane Surface(92) ={92}; -Curve Loop(93) = { -106, 228, 107, -218}; -Plane Surface(93) ={93}; -Curve Loop(94) = { -107, -229, 108, 219}; -Plane Surface(94) ={94}; -Curve Loop(95) = { -108, 230, 109, -220}; -Plane Surface(95) ={95}; -Curve Loop(96) = { -109, -231, 110, 221}; -Plane Surface(96) ={96}; -Curve Loop(97) = { -110, 232, 111, -222}; -Plane Surface(97) ={97}; -Curve Loop(98) = { -111, -233, 112, 223}; -Plane Surface(98) ={98}; -Curve Loop(99) = { -112, 234, 113, -224}; -Plane Surface(99) ={99}; -Curve Loop(100) = { -113, -235, 114, 225}; -Plane Surface(100) ={100}; - - -// "Outside area" -//Curve Loop(101) = {1, 2, 3, 4}; -//+ -//Curve Loop(102) = {5, -16, 27, -38, 49, -60, 71, -82, 93, -104, 226, 227, 228, 229, 230, 231, 232, 233, 234, 235, 114, -103, 92, -81, 70, -59, 48, -37, 26, -15, -135, -134, -133, -132, -131, -130, -129, -128, -127, -126}; -//+ -//Plane Surface(101) = {101, 102}; - -//+ -Transfinite Surface {1}; -Transfinite Surface {2}; -Transfinite Surface {3}; -Transfinite Surface {4}; -Transfinite Surface {5}; -Transfinite Surface {6}; -Transfinite Surface {7}; -Transfinite Surface {8}; -Transfinite Surface {9}; -Transfinite Surface {10}; -Transfinite Surface {11}; -Transfinite Surface {12}; -Transfinite Surface {13}; -Transfinite Surface {14}; -Transfinite Surface {15}; -Transfinite Surface {16}; -Transfinite Surface {17}; -Transfinite Surface {18}; -Transfinite Surface {19}; -Transfinite Surface {20}; -Transfinite Surface {21}; -Transfinite Surface {22}; -Transfinite Surface {23}; -Transfinite Surface {24}; -Transfinite Surface {25}; -Transfinite Surface {26}; -Transfinite Surface {27}; -Transfinite Surface {28}; -Transfinite Surface {29}; -Transfinite Surface {30}; -Transfinite Surface {31}; -Transfinite Surface {32}; -Transfinite Surface {33}; -Transfinite Surface {34}; -Transfinite Surface {35}; -Transfinite Surface {36}; -Transfinite Surface {37}; -Transfinite Surface {38}; -Transfinite Surface {39}; -Transfinite Surface {40}; -Transfinite Surface {41}; -Transfinite Surface {42}; -Transfinite Surface {43}; -Transfinite Surface {44}; -Transfinite Surface {45}; -Transfinite Surface {46}; -Transfinite Surface {47}; -Transfinite Surface {48}; -Transfinite Surface {49}; -Transfinite Surface {50}; -Transfinite Surface {51}; -Transfinite Surface {52}; -Transfinite Surface {53}; -Transfinite Surface {54}; -Transfinite Surface {55}; -Transfinite Surface {56}; -Transfinite Surface {57}; -Transfinite Surface {58}; -Transfinite Surface {59}; -Transfinite Surface {60}; -Transfinite Surface {61}; -Transfinite Surface {62}; -Transfinite Surface {63}; -Transfinite Surface {64}; -Transfinite Surface {65}; -Transfinite Surface {66}; -Transfinite Surface {67}; -Transfinite Surface {68}; -Transfinite Surface {69}; -Transfinite Surface {70}; -Transfinite Surface {71}; -Transfinite Surface {72}; -Transfinite Surface {73}; -Transfinite Surface {74}; -Transfinite Surface {75}; -Transfinite Surface {76}; -Transfinite Surface {77}; -Transfinite Surface {78}; -Transfinite Surface {79}; -Transfinite Surface {80}; -Transfinite Surface {81}; -Transfinite Surface {82}; -Transfinite Surface {83}; -Transfinite Surface {84}; -Transfinite Surface {85}; -Transfinite Surface {86}; -Transfinite Surface {87}; -Transfinite Surface {88}; -Transfinite Surface {89}; -Transfinite Surface {90}; -Transfinite Surface {91}; -Transfinite Surface {92}; -Transfinite Surface {93}; -Transfinite Surface {94}; -Transfinite Surface {95}; -Transfinite Surface {96}; -Transfinite Surface {97}; -Transfinite Surface {98}; -Transfinite Surface {99}; -Transfinite Surface {100}; -//+ -Transfinite Curve {1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103} = n_recombine Using Progression n_prog; -// add outer perimeter with half as many cells and no Progression -Transfinite Curve {104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15} = n_recombine / 4 Using Progression 1; -//+ -Transfinite Curve { 127, 137, 147, 157, 167, 177, 187, 197, 207, 217, 227, 128, 138, 148, 158, 168, 178, 188, 198, 208, 218, 228, 129, 139, 149, 159, 169, 179, 189, 199, 209, 219, 229, 130, 140, 150, 160, 170, 180, 190, 200, 210, 220, 230, 131, 141, 151, 161, 171, 181, 191, 201, 211, 221, 231, 132, 142, 152, 162, 172, 182, 192, 202, 212, 222, 232, 133, 143, 153, 163, 173, 183, 193, 203, 213, 223, 233, 134, 144, 154, 164, 174, 184, 194, 204, 214, 224, 234} = n_recombine Using Progression n_prog; -// add outer perimeter with half as many cells and no Progression -Transfinite Curve {126, 136, 146, 156, 166, 176, 186, 196, 206, 216, 226, 135, 145, 155, 165, 175, 185, 195, 205, 215, 225, 235} = n_recombine / 4 Using Progression 1; -//+ -Line(236) = {1, 126}; -//+ -Line(237) = {126, 127}; -//+ -Line(238) = {127, 128}; -//+ -Line(239) = {128, 129}; -//+ -Line(240) = {129, 130}; -//+ -Line(241) = {130, 131}; -//+ -Line(242) = {131, 132}; -//+ -Line(243) = {132, 133}; -//+ -Line(244) = {133, 134}; -//+ -Line(245) = {134, 135}; -//+ -Line(246) = {135, 136}; -//+ -Line(247) = {136, 2}; -//+ -Line(248) = {2, 159}; -//+ -Line(249) = {159, 160}; -//+ -Line(250) = {160, 161}; -//+ -Line(251) = {161, 162}; -//+ -Line(252) = {162, 163}; -//+ -Line(253) = {163, 164}; -//+ -Line(254) = {164, 165}; -//+ -Line(255) = {165, 166}; -//+ -Line(256) = {166, 167}; -//+ -Line(257) = {167, 168}; -//+ -Line(258) = {168, 169}; -//+ -Line(259) = {169, 4}; -//+ -Line(260) = {4, 147}; -//+ -Line(261) = {147, 146}; -//+ -Line(262) = {146, 145}; -//+ -Line(263) = {145, 144}; -//+ -Line(264) = {144, 143}; -//+ -Line(265) = {143, 142}; -//+ -Line(266) = {142, 141}; -//+ -Line(267) = {141, 140}; -//+ -Line(268) = {140, 139}; -//+ -Line(269) = {139, 138}; -//+ -Line(270) = {138, 137}; -//+ -Line(271) = {137, 3}; -//+ -Line(272) = {3, 158}; -//+ -Line(273) = {158, 157}; -//+ -Line(274) = {157, 156}; -//+ -Line(275) = {156, 155}; -//+ -Line(276) = {155, 154}; -//+ -Line(277) = {154, 153}; -//+ -Line(278) = {153, 152}; -//+ -Line(279) = {152, 151}; -//+ -Line(280) = {151, 150}; -//+ -Line(281) = {150, 149}; -//+ -Line(282) = {149, 148}; -//+ -Line(283) = {148, 1}; -//+ -Line(284) = {126, 5}; -//+ -Line(285) = {127, 16}; -//+ -Line(286) = {128, 27}; -//+ -Line(287) = {129, 38}; -//+ -Line(288) = {130, 49}; -//+ -Line(289) = {131, 60}; -//+ -Line(290) = {132, 71}; -//+ -Line(291) = {133, 82}; -//+ -Line(292) = {134, 93}; -//+ -Line(293) = {135, 104}; -//+ -Line(294) = {136, 115}; -//+ -Line(295) = {159, 115}; -//+ -Line(296) = {160, 116}; -//+ -Line(297) = {161, 117}; -//+ -Line(298) = {162, 118}; -//+ -Line(299) = {163, 119}; -//+ -Line(300) = {164, 120}; -//+ -Line(301) = {165, 121}; -//+ -Line(302) = {166, 122}; -//+ -Line(303) = {167, 123}; -//+ -Line(304) = {168, 124}; -//+ -Line(305) = {169, 125}; -//+ -Line(306) = {147, 125}; -//+ -Line(307) = {146, 114}; -//+ -Line(308) = {145, 103}; -//+ -Line(309) = {144, 92}; -//+ -Line(310) = {143, 81}; -//+ -Line(311) = {142, 70}; -//+ -Line(312) = {141, 59}; -//+ -Line(313) = {140, 48}; -//+ -Line(314) = {139, 37}; -//+ -Line(315) = {138, 26}; -//+ -Line(316) = {137, 15}; -//+ -Line(317) = {158, 15}; -//+ -Line(318) = {157, 14}; -//+ -Line(319) = {156, 13}; -//+ -Line(320) = {155, 12}; -//+ -Line(321) = {154, 11}; -//+ -Line(322) = {153, 10}; -//+ -Line(323) = {152, 9}; -//+ -Line(324) = {151, 8}; -//+ -Line(325) = {150, 7}; -//+ -Line(326) = {149, 6}; -//+ -Line(327) = {148, 5}; -//+ -Curve Loop(103) = {236, 284, -327, 283}; -//+ -Plane Surface(102) = {103}; -//+ -Curve Loop(104) = {237, 285, -5, -284}; -//+ -Plane Surface(103) = {104}; -//+ -Curve Loop(105) = {238, 286, 16, -285}; -//+ -Plane Surface(104) = {105}; -//+ -Curve Loop(106) = {239, 287, -27, -286}; -//+ -Plane Surface(105) = {106}; -//+ -Curve Loop(107) = {240, 288, 38, -287}; -//+ -Plane Surface(106) = {107}; -//+ -Curve Loop(108) = {241, 289, -49, -288}; -//+ -Plane Surface(107) = {108}; -//+ -Curve Loop(109) = {242, 290, 60, -289}; -//+ -Plane Surface(108) = {109}; -//+ -Curve Loop(110) = {243, 291, -71, -290}; -//+ -Plane Surface(109) = {110}; -//+ -Curve Loop(111) = {244, 292, 82, -291}; -//+ -Plane Surface(110) = {111}; -//+ -Curve Loop(112) = {245, 293, -93, -292}; -//+ -Plane Surface(111) = {112}; -//+ -Curve Loop(113) = {246, 294, 104, -293}; -//+ -Plane Surface(112) = {113}; -//+ -Curve Loop(114) = {247, 248, 295, -294}; -//+ -Plane Surface(113) = {114}; -//+ -Curve Loop(115) = {249, 296, -226, -295}; -//+ -Plane Surface(114) = {115}; -//+ -Curve Loop(116) = {250, 297, 227, -296}; -//+ -Plane Surface(115) = {116}; -//+ -Curve Loop(117) = {251, 298, -228, -297}; -//+ -Plane Surface(116) = {117}; -//+ -Curve Loop(118) = {252, 299, 229, -298}; -//+ -Plane Surface(117) = {118}; -//+ -Curve Loop(119) = {253, 300, -230, -299}; -//+ -Plane Surface(118) = {119}; -//+ -Curve Loop(120) = {254, 301, 231, -300}; -//+ -Plane Surface(119) = {120}; -//+ -Curve Loop(121) = {255, 302, -232, -301}; -//+ -Plane Surface(120) = {121}; -//+ -Curve Loop(122) = {256, 303, 233, -302}; -//+ -Plane Surface(121) = {122}; -//+ -Curve Loop(123) = {257, 304, -234, -303}; -//+ -Plane Surface(122) = {123}; -//+ -Curve Loop(124) = {258, 305, 235, -304}; -//+ -Plane Surface(123) = {124}; -//+ -Curve Loop(125) = {259, 260, 306, -305}; -//+ -Plane Surface(124) = {125}; -//+ -Curve Loop(126) = {261, 307, -114, -306}; -//+ -Plane Surface(125) = {126}; -//+ -Curve Loop(127) = {262, 308, 103, -307}; -//+ -Plane Surface(126) = {127}; -//+ -Curve Loop(128) = {263, 309, -92, -308}; -//+ -Plane Surface(127) = {128}; -//+ -Curve Loop(129) = {264, 310, 81, -309}; -//+ -Plane Surface(128) = {129}; -//+ -Curve Loop(130) = {265, 311, -70, -310}; -//+ -Plane Surface(129) = {130}; -//+ -Curve Loop(131) = {266, 312, 59, -311}; -//+ -Plane Surface(130) = {131}; -//+ -Curve Loop(132) = {267, 313, -48, -312}; -//+ -Plane Surface(131) = {132}; -//+ -Curve Loop(133) = {268, 314, 37, -313}; -//+ -Plane Surface(132) = {133}; -//+ -Curve Loop(134) = {269, 315, -26, -314}; -//+ -Plane Surface(133) = {134}; -//+ -Curve Loop(135) = {270, 316, 15, -315}; -//+ -Plane Surface(134) = {135}; -//+ -Curve Loop(136) = {271, 272, 317, -316}; -//+ -Plane Surface(135) = {136}; -//+ -Curve Loop(137) = {273, 318, -135, -317}; -//+ -Plane Surface(136) = {137}; -//+ -Curve Loop(138) = {274, 319, 134, -318}; -//+ -Plane Surface(137) = {138}; -//+ -Curve Loop(139) = {275, 320, -133, -319}; -//+ -Plane Surface(138) = {139}; -//+ -Curve Loop(140) = {276, 321, 132, -320}; -//+ -Plane Surface(139) = {140}; -//+ -Curve Loop(141) = {277, 322, -131, -321}; -//+ -Plane Surface(140) = {141}; -//+ -Curve Loop(142) = {278, 323, 130, -322}; -//+ -Plane Surface(141) = {142}; -//+ -Curve Loop(143) = {279, 324, -129, -323}; -//+ -Plane Surface(142) = {143}; -//+ -Curve Loop(144) = {280, 325, 128, -324}; -//+ -Plane Surface(143) = {144}; -//+ -Curve Loop(145) = {281, 326, -127, -325}; -//+ -Plane Surface(144) = {145}; -//+ -Curve Loop(146) = {282, 327, 126, -326}; -//+ -Plane Surface(145) = {146}; -//+ -//+ -Transfinite Surface {102}; -//+ -Transfinite Surface {103}; -//+ -Transfinite Surface {104}; -//+ -Transfinite Surface {105}; -//+ -Transfinite Surface {106}; -//+ -Transfinite Surface {107}; -//+ -Transfinite Surface {108}; -//+ -Transfinite Surface {109}; -//+ -Transfinite Surface {110}; -//+ -Transfinite Surface {111}; -//+ -Transfinite Surface {112}; -//+ -Transfinite Surface {113}; -//+ -Transfinite Surface {114}; -//+ -Transfinite Surface {115}; -//+ -Transfinite Surface {116}; -//+ -Transfinite Surface {117}; -//+ -Transfinite Surface {118}; -//+ -Transfinite Surface {119}; -//+ -Transfinite Surface {120}; -//+ -Transfinite Surface {121}; -//+ -Transfinite Surface {122}; -//+ -Transfinite Surface {123}; -//+ -Transfinite Surface {124}; -//+ -Transfinite Surface {125}; -//+ -Transfinite Surface {126}; -//+ -Transfinite Surface {127}; -//+ -Transfinite Surface {128}; -//+ -Transfinite Surface {129}; -//+ -Transfinite Surface {130}; -//+ -Transfinite Surface {131}; -//+ -Transfinite Surface {132}; -//+ -Transfinite Surface {133}; -//+ -Transfinite Surface {134}; -//+ -Transfinite Surface {135}; -//+ -Transfinite Surface {136}; -//+ -Transfinite Surface {137}; -//+ -Transfinite Surface {138}; -//+ -Transfinite Surface {139}; -//+ -Transfinite Surface {140}; -//+ -Transfinite Surface {141}; -//+ -Transfinite Surface {142}; -//+ -Transfinite Surface {143}; -//+ -Transfinite Surface {144}; -//+ -Transfinite Surface {145}; -//+ -Transfinite Curve { 238, 239, 240, 241, 242, 243, 244, 245, 250, 251, 252, 253, 254, 255, 256, 257, 262, 263, 264, 265, 266, 267, 268, 269, 274, 275, 276, 277, 278, 279, 280, 281} = n_recombine Using Progression 1; -Transfinite Curve {237,246,249, 258,261,270,273,282} = n_recombine / 4 Using Progression 1; - - -Physical Curve("void",328) = {236, 237, 238, 239, 240, 241, 242, 243, 244, 245, 246, 247, 248, 249, 250, 251, 252, 253, 254, 255, 256, 257, 258, 259, 260, 261, 262, 263, 264, 265, 266, 267, 268, 269, 270, 271, 272, 273, 274, 275, 276, 277, 278, 279, 280, 281, 282, 283}; - -Recombine Surface "*"; -// Define meshing options -//Mesh.Algorithm = 6; // Specify meshing algorithm (e.g., Delaunay) -//Mesh.ElementOrder = 1; // Specify element order -// Generate the mesh -//Mesh 2; \ No newline at end of file diff --git a/examples_hpc/lattice/mesh/lattice_rectangular_n40.geo b/examples_hpc/lattice/mesh/lattice_rectangular_n40.geo deleted file mode 100644 index ccc4e41f..00000000 --- a/examples_hpc/lattice/mesh/lattice_rectangular_n40.geo +++ /dev/null @@ -1,1193 +0,0 @@ -cl_fine = 0.1; -n_recombine = 40; -n_prog = 1.1; -Point(1) = {-3.5, -3.5, 0, cl_fine}; -Point(2) = {3.5, -3.5, 0, cl_fine}; -Point(3) = {-3.5, 3.5, 0, cl_fine}; -Point(4) = {3.5, 3.5, 0, cl_fine}; - -// Inner grid -Point(5) = { -2.5, -2.5, 0, cl_fine }; -Point(6) = { -2.5, -2.0, 0, cl_fine }; -Point(7) = { -2.5, -1.5, 0, cl_fine }; -Point(8) = { -2.5, -1.0, 0, cl_fine }; -Point(9) = { -2.5, -0.5, 0, cl_fine }; -Point(10) = { -2.5, 0.0, 0, cl_fine }; -Point(11) = { -2.5, 0.5, 0, cl_fine }; -Point(12) = { -2.5, 1.0, 0, cl_fine }; -Point(13) = { -2.5, 1.5, 0, cl_fine }; -Point(14) = { -2.5, 2.0, 0, cl_fine }; -Point(15) = { -2.5, 2.5, 0, cl_fine }; -Point(16) = { -2.0, -2.5, 0, cl_fine }; -Point(17) = { -2.0, -2.0, 0, cl_fine }; -Point(18) = { -2.0, -1.5, 0, cl_fine }; -Point(19) = { -2.0, -1.0, 0, cl_fine }; -Point(20) = { -2.0, -0.5, 0, cl_fine }; -Point(21) = { -2.0, 0.0, 0, cl_fine }; -Point(22) = { -2.0, 0.5, 0, cl_fine }; -Point(23) = { -2.0, 1.0, 0, cl_fine }; -Point(24) = { -2.0, 1.5, 0, cl_fine }; -Point(25) = { -2.0, 2.0, 0, cl_fine }; -Point(26) = { -2.0, 2.5, 0, cl_fine }; -Point(27) = { -1.5, -2.5, 0, cl_fine }; -Point(28) = { -1.5, -2.0, 0, cl_fine }; -Point(29) = { -1.5, -1.5, 0, cl_fine }; -Point(30) = { -1.5, -1.0, 0, cl_fine }; -Point(31) = { -1.5, -0.5, 0, cl_fine }; -Point(32) = { -1.5, 0.0, 0, cl_fine }; -Point(33) = { -1.5, 0.5, 0, cl_fine }; -Point(34) = { -1.5, 1.0, 0, cl_fine }; -Point(35) = { -1.5, 1.5, 0, cl_fine }; -Point(36) = { -1.5, 2.0, 0, cl_fine }; -Point(37) = { -1.5, 2.5, 0, cl_fine }; -Point(38) = { -1.0, -2.5, 0, cl_fine }; -Point(39) = { -1.0, -2.0, 0, cl_fine }; -Point(40) = { -1.0, -1.5, 0, cl_fine }; -Point(41) = { -1.0, -1.0, 0, cl_fine }; -Point(42) = { -1.0, -0.5, 0, cl_fine }; -Point(43) = { -1.0, 0.0, 0, cl_fine }; -Point(44) = { -1.0, 0.5, 0, cl_fine }; -Point(45) = { -1.0, 1.0, 0, cl_fine }; -Point(46) = { -1.0, 1.5, 0, cl_fine }; -Point(47) = { -1.0, 2.0, 0, cl_fine }; -Point(48) = { -1.0, 2.5, 0, cl_fine }; -Point(49) = { -0.5, -2.5, 0, cl_fine }; -Point(50) = { -0.5, -2.0, 0, cl_fine }; -Point(51) = { -0.5, -1.5, 0, cl_fine }; -Point(52) = { -0.5, -1.0, 0, cl_fine }; -Point(53) = { -0.5, -0.5, 0, cl_fine }; -Point(54) = { -0.5, 0.0, 0, cl_fine }; -Point(55) = { -0.5, 0.5, 0, cl_fine }; -Point(56) = { -0.5, 1.0, 0, cl_fine }; -Point(57) = { -0.5, 1.5, 0, cl_fine }; -Point(58) = { -0.5, 2.0, 0, cl_fine }; -Point(59) = { -0.5, 2.5, 0, cl_fine }; -Point(60) = { 0.0, -2.5, 0, cl_fine }; -Point(61) = { 0.0, -2.0, 0, cl_fine }; -Point(62) = { 0.0, -1.5, 0, cl_fine }; -Point(63) = { 0.0, -1.0, 0, cl_fine }; -Point(64) = { 0.0, -0.5, 0, cl_fine }; -Point(65) = { 0.0, 0.0, 0, cl_fine }; -Point(66) = { 0.0, 0.5, 0, cl_fine }; -Point(67) = { 0.0, 1.0, 0, cl_fine }; -Point(68) = { 0.0, 1.5, 0, cl_fine }; -Point(69) = { 0.0, 2.0, 0, cl_fine }; -Point(70) = { 0.0, 2.5, 0, cl_fine }; -Point(71) = { 0.5, -2.5, 0, cl_fine }; -Point(72) = { 0.5, -2.0, 0, cl_fine }; -Point(73) = { 0.5, -1.5, 0, cl_fine }; -Point(74) = { 0.5, -1.0, 0, cl_fine }; -Point(75) = { 0.5, -0.5, 0, cl_fine }; -Point(76) = { 0.5, 0.0, 0, cl_fine }; -Point(77) = { 0.5, 0.5, 0, cl_fine }; -Point(78) = { 0.5, 1.0, 0, cl_fine }; -Point(79) = { 0.5, 1.5, 0, cl_fine }; -Point(80) = { 0.5, 2.0, 0, cl_fine }; -Point(81) = { 0.5, 2.5, 0, cl_fine }; -Point(82) = { 1.0, -2.5, 0, cl_fine }; -Point(83) = { 1.0, -2.0, 0, cl_fine }; -Point(84) = { 1.0, -1.5, 0, cl_fine }; -Point(85) = { 1.0, -1.0, 0, cl_fine }; -Point(86) = { 1.0, -0.5, 0, cl_fine }; -Point(87) = { 1.0, 0.0, 0, cl_fine }; -Point(88) = { 1.0, 0.5, 0, cl_fine }; -Point(89) = { 1.0, 1.0, 0, cl_fine }; -Point(90) = { 1.0, 1.5, 0, cl_fine }; -Point(91) = { 1.0, 2.0, 0, cl_fine }; -Point(92) = { 1.0, 2.5, 0, cl_fine }; -Point(93) = { 1.5, -2.5, 0, cl_fine }; -Point(94) = { 1.5, -2.0, 0, cl_fine }; -Point(95) = { 1.5, -1.5, 0, cl_fine }; -Point(96) = { 1.5, -1.0, 0, cl_fine }; -Point(97) = { 1.5, -0.5, 0, cl_fine }; -Point(98) = { 1.5, 0.0, 0, cl_fine }; -Point(99) = { 1.5, 0.5, 0, cl_fine }; -Point(100) = { 1.5, 1.0, 0, cl_fine }; -Point(101) = { 1.5, 1.5, 0, cl_fine }; -Point(102) = { 1.5, 2.0, 0, cl_fine }; -Point(103) = { 1.5, 2.5, 0, cl_fine }; -Point(104) = { 2.0, -2.5, 0, cl_fine }; -Point(105) = { 2.0, -2.0, 0, cl_fine }; -Point(106) = { 2.0, -1.5, 0, cl_fine }; -Point(107) = { 2.0, -1.0, 0, cl_fine }; -Point(108) = { 2.0, -0.5, 0, cl_fine }; -Point(109) = { 2.0, 0.0, 0, cl_fine }; -Point(110) = { 2.0, 0.5, 0, cl_fine }; -Point(111) = { 2.0, 1.0, 0, cl_fine }; -Point(112) = { 2.0, 1.5, 0, cl_fine }; -Point(113) = { 2.0, 2.0, 0, cl_fine }; -Point(114) = { 2.0, 2.5, 0, cl_fine }; -Point(115) = { 2.5, -2.5, 0, cl_fine }; -Point(116) = { 2.5, -2.0, 0, cl_fine }; -Point(117) = { 2.5, -1.5, 0, cl_fine }; -Point(118) = { 2.5, -1.0, 0, cl_fine }; -Point(119) = { 2.5, -0.5, 0, cl_fine }; -Point(120) = { 2.5, 0.0, 0, cl_fine }; -Point(121) = { 2.5, 0.5, 0, cl_fine }; -Point(122) = { 2.5, 1.0, 0, cl_fine }; -Point(123) = { 2.5, 1.5, 0, cl_fine }; -Point(124) = { 2.5, 2.0, 0, cl_fine }; -Point(125) = { 2.5, 2.5, 0, cl_fine }; - -// helper boundary points -Point(126) = { -2.5, -3.5, 0, cl_fine }; -Point(127) = { -2.0, -3.5, 0, cl_fine }; -Point(128) = { -1.5, -3.5, 0, cl_fine }; -Point(129) = { -1.0, -3.5, 0, cl_fine }; -Point(130) = { -0.5, -3.5, 0, cl_fine }; -Point(131) = { 0.0, -3.5, 0, cl_fine }; -Point(132) = { 0.5, -3.5, 0, cl_fine }; -Point(133) = { 1.0, -3.5, 0, cl_fine }; -Point(134) = { 1.5, -3.5, 0, cl_fine }; -Point(135) = { 2.0, -3.5, 0, cl_fine }; -Point(136) = { 2.5, -3.5, 0, cl_fine }; -Point(137) = { -2.5, 3.5, 0, cl_fine }; -Point(138) = { -2.0, 3.5, 0, cl_fine }; -Point(139) = { -1.5, 3.5, 0, cl_fine }; -Point(140) = { -1.0, 3.5, 0, cl_fine }; -Point(141) = { -0.5, 3.5, 0, cl_fine }; -Point(142) = { 0.0, 3.5, 0, cl_fine }; -Point(143) = { 0.5, 3.5, 0, cl_fine }; -Point(144) = { 1.0, 3.5, 0, cl_fine }; -Point(145) = { 1.5, 3.5, 0, cl_fine }; -Point(146) = { 2.0, 3.5, 0, cl_fine }; -Point(147) = { 2.5, 3.5, 0, cl_fine }; -Point(148) = { -3.5, -2.5, 0, cl_fine }; -Point(149) = { -3.5, -2.0, 0, cl_fine }; -Point(150) = { -3.5, -1.5, 0, cl_fine }; -Point(151) = { -3.5, -1.0, 0, cl_fine }; -Point(152) = { -3.5, -0.5, 0, cl_fine }; -Point(153) = { -3.5, 0.0, 0, cl_fine }; -Point(154) = { -3.5, 0.5, 0, cl_fine }; -Point(155) = { -3.5, 1.0, 0, cl_fine }; -Point(156) = { -3.5, 1.5, 0, cl_fine }; -Point(157) = { -3.5, 2.0, 0, cl_fine }; -Point(158) = { -3.5, 2.5, 0, cl_fine }; -Point(159) = { 3.5, -2.5, 0, cl_fine }; -Point(160) = { 3.5, -2.0, 0, cl_fine }; -Point(161) = { 3.5, -1.5, 0, cl_fine }; -Point(162) = { 3.5, -1.0, 0, cl_fine }; -Point(163) = { 3.5, -0.5, 0, cl_fine }; -Point(164) = { 3.5, 0.0, 0, cl_fine }; -Point(165) = { 3.5, 0.5, 0, cl_fine }; -Point(166) = { 3.5, 1.0, 0, cl_fine }; -Point(167) = { 3.5, 1.5, 0, cl_fine }; -Point(168) = { 3.5, 2.0, 0, cl_fine }; -Point(169) = { 3.5, 2.5, 0, cl_fine }; - -// Horizontal and vertical lines in inner grid -Line(5) = { 5, 16}; -Line(6) = { 6, 17}; -Line(7) = { 7, 18}; -Line(8) = { 8, 19}; -Line(9) = { 9, 20}; -Line(10) = { 10, 21}; -Line(11) = { 11, 22}; -Line(12) = { 12, 23}; -Line(13) = { 13, 24}; -Line(14) = { 14, 25}; -Line(15) = { 15, 26}; -Line(16) = { 27, 16}; -Line(17) = { 28, 17}; -Line(18) = { 29, 18}; -Line(19) = { 30, 19}; -Line(20) = { 31, 20}; -Line(21) = { 32, 21}; -Line(22) = { 33, 22}; -Line(23) = { 34, 23}; -Line(24) = { 35, 24}; -Line(25) = { 36, 25}; -Line(26) = { 37, 26}; -Line(27) = { 27, 38}; -Line(28) = { 28, 39}; -Line(29) = { 29, 40}; -Line(30) = { 30, 41}; -Line(31) = { 31, 42}; -Line(32) = { 32, 43}; -Line(33) = { 33, 44}; -Line(34) = { 34, 45}; -Line(35) = { 35, 46}; -Line(36) = { 36, 47}; -Line(37) = { 37, 48}; -Line(38) = { 49, 38}; -Line(39) = { 50, 39}; -Line(40) = { 51, 40}; -Line(41) = { 52, 41}; -Line(42) = { 53, 42}; -Line(43) = { 54, 43}; -Line(44) = { 55, 44}; -Line(45) = { 56, 45}; -Line(46) = { 57, 46}; -Line(47) = { 58, 47}; -Line(48) = { 59, 48}; -Line(49) = { 49, 60}; -Line(50) = { 50, 61}; -Line(51) = { 51, 62}; -Line(52) = { 52, 63}; -Line(53) = { 53, 64}; -Line(54) = { 54, 65}; -Line(55) = { 55, 66}; -Line(56) = { 56, 67}; -Line(57) = { 57, 68}; -Line(58) = { 58, 69}; -Line(59) = { 59, 70}; -Line(60) = { 71, 60}; -Line(61) = { 72, 61}; -Line(62) = { 73, 62}; -Line(63) = { 74, 63}; -Line(64) = { 75, 64}; -Line(65) = { 76, 65}; -Line(66) = { 77, 66}; -Line(67) = { 78, 67}; -Line(68) = { 79, 68}; -Line(69) = { 80, 69}; -Line(70) = { 81, 70}; -Line(71) = { 71, 82}; -Line(72) = { 72, 83}; -Line(73) = { 73, 84}; -Line(74) = { 74, 85}; -Line(75) = { 75, 86}; -Line(76) = { 76, 87}; -Line(77) = { 77, 88}; -Line(78) = { 78, 89}; -Line(79) = { 79, 90}; -Line(80) = { 80, 91}; -Line(81) = { 81, 92}; -Line(82) = { 93, 82}; -Line(83) = { 94, 83}; -Line(84) = { 95, 84}; -Line(85) = { 96, 85}; -Line(86) = { 97, 86}; -Line(87) = { 98, 87}; -Line(88) = { 99, 88}; -Line(89) = { 100, 89}; -Line(90) = { 101, 90}; -Line(91) = { 102, 91}; -Line(92) = { 103, 92}; -Line(93) = { 93, 104}; -Line(94) = { 94, 105}; -Line(95) = { 95, 106}; -Line(96) = { 96, 107}; -Line(97) = { 97, 108}; -Line(98) = { 98, 109}; -Line(99) = { 99, 110}; -Line(100) = { 100, 111}; -Line(101) = { 101, 112}; -Line(102) = { 102, 113}; -Line(103) = { 103, 114}; -Line(104) = { 115, 104}; -Line(105) = { 116, 105}; -Line(106) = { 117, 106}; -Line(107) = { 118, 107}; -Line(108) = { 119, 108}; -Line(109) = { 120, 109}; -Line(110) = { 121, 110}; -Line(111) = { 122, 111}; -Line(112) = { 123, 112}; -Line(113) = { 124, 113}; -Line(114) = { 125, 114}; -Line(126) = { 5, 6}; -Line(127) = { 7, 6}; -Line(128) = { 7, 8}; -Line(129) = { 9, 8}; -Line(130) = { 9, 10}; -Line(131) = { 11, 10}; -Line(132) = { 11, 12}; -Line(133) = { 13, 12}; -Line(134) = { 13, 14}; -Line(135) = { 15, 14}; -Line(136) = { 16, 17}; -Line(137) = { 18, 17}; -Line(138) = { 18, 19}; -Line(139) = { 20, 19}; -Line(140) = { 20, 21}; -Line(141) = { 22, 21}; -Line(142) = { 22, 23}; -Line(143) = { 24, 23}; -Line(144) = { 24, 25}; -Line(145) = { 26, 25}; - -Line(146) = { 27, 28}; -Line(147) = { 29, 28}; -Line(148) = { 29, 30}; -Line(149) = { 31, 30}; -Line(150) = { 31, 32}; -Line(151) = { 33, 32}; -Line(152) = { 33, 34}; -Line(153) = { 35, 34}; -Line(154) = { 35, 36}; -Line(155) = { 37, 36}; - -Line(156) = { 38, 39}; -Line(157) = { 40, 39}; -Line(158) = { 40, 41}; -Line(159) = { 42, 41}; -Line(160) = { 42, 43}; -Line(161) = { 44, 43}; -Line(162) = { 44, 45}; -Line(163) = { 46, 45}; -Line(164) = { 46, 47}; -Line(165) = { 48, 47}; - -Line(166) = { 49, 50}; -Line(167) = { 51, 50}; -Line(168) = { 51, 52}; -Line(169) = { 53, 52}; -Line(170) = { 53, 54}; -Line(171) = { 55, 54}; -Line(172) = { 55, 56}; -Line(173) = { 57, 56}; -Line(174) = { 57, 58}; -Line(175) = { 59, 58}; - -Line(176) = { 60, 61}; -Line(177) = { 62, 61}; -Line(178) = { 62, 63}; -Line(179) = { 64, 63}; -Line(180) = { 64, 65}; -Line(181) = { 66, 65}; -Line(182) = { 66, 67}; -Line(183) = { 68, 67}; -Line(184) = { 68, 69}; -Line(185) = { 70, 69}; - -Line(186) = { 71, 72}; -Line(187) = { 73, 72}; -Line(188) = { 73, 74}; -Line(189) = { 75, 74}; -Line(190) = { 75, 76}; -Line(191) = { 77, 76}; -Line(192) = { 77, 78}; -Line(193) = { 79, 78}; -Line(194) = { 79, 80}; -Line(195) = { 81, 80}; - -Line(196) = { 82, 83}; -Line(197) = { 84, 83}; -Line(198) = { 84, 85}; -Line(199) = { 86, 85}; -Line(200) = { 86, 87}; -Line(201) = { 88, 87}; -Line(202) = { 88, 89}; -Line(203) = { 90, 89}; -Line(204) = { 90, 91}; -Line(205) = { 92, 91}; - -Line(206) = { 93, 94}; -Line(207) = { 95, 94}; -Line(208) = { 95, 96}; -Line(209) = { 97, 96}; -Line(210) = { 97, 98}; -Line(211) = { 99, 98}; -Line(212) = { 99, 100}; -Line(213) = { 101, 100}; -Line(214) = { 101, 102}; -Line(215) = { 103, 102}; - - -Line(216) = { 104, 105}; -Line(217) = { 106, 105}; -Line(218) = { 106, 107}; -Line(219) = { 108, 107}; -Line(220) = { 108, 109}; -Line(221) = { 110, 109}; -Line(222) = { 110, 111}; -Line(223) = { 112, 111}; -Line(224) = { 112, 113}; -Line(225) = { 114, 113}; - - -Line(226) = { 115, 116}; -Line(227) = { 117, 116}; -Line(228) = { 117, 118}; -Line(229) = { 119, 118}; -Line(230) = { 119, 120}; -Line(231) = { 121, 120}; -Line(232) = { 121, 122}; -Line(233) = { 123, 122}; -Line(234) = { 123, 124}; -Line(235) = { 125, 124}; - - -// Define inner surfaces -//+ -Curve Loop(1) = { 5, 136, -6, -126}; -Plane Surface(1) ={1}; -Curve Loop(2) = { 6, -137, -7, 127}; -Plane Surface(2) ={2}; -Curve Loop(3) = { 7, 138, -8, -128}; -Plane Surface(3) ={3}; -Curve Loop(4) = { 8, -139, -9, 129}; -Plane Surface(4) ={4}; -Curve Loop(5) = { 9, 140, -10, -130}; -Plane Surface(5) ={5}; -Curve Loop(6) = { 10, -141, -11, 131}; -Plane Surface(6) ={6}; -Curve Loop(7) = { 11, 142, -12, -132}; -Plane Surface(7) ={7}; -Curve Loop(8) = { 12, -143, -13, 133}; -Plane Surface(8) ={8}; -Curve Loop(9) = { 13, 144, -14, -134}; -Plane Surface(9) ={9}; -Curve Loop(10) = { 14, -145, -15, 135}; -Plane Surface(10) ={10}; -Curve Loop(11) = { -16, 146, 17, -136}; -Plane Surface(11) ={11}; -Curve Loop(12) = { -17, -147, 18, 137}; -Plane Surface(12) ={12}; -Curve Loop(13) = { -18, 148, 19, -138}; -Plane Surface(13) ={13}; -Curve Loop(14) = { -19, -149, 20, 139}; -Plane Surface(14) ={14}; -Curve Loop(15) = { -20, 150, 21, -140}; -Plane Surface(15) ={15}; -Curve Loop(16) = { -21, -151, 22, 141}; -Plane Surface(16) ={16}; -Curve Loop(17) = { -22, 152, 23, -142}; -Plane Surface(17) ={17}; -Curve Loop(18) = { -23, -153, 24, 143}; -Plane Surface(18) ={18}; -Curve Loop(19) = {-24, 154, 25, -144}; -Plane Surface(19) ={19}; -Curve Loop(20) = { -25, -155, 26, 145}; -Plane Surface(20) ={20}; -Curve Loop(21) = { 27, 156, -28, -146}; -Plane Surface(21) ={21}; -Curve Loop(22) = { 28, -157, -29, 147}; -Plane Surface(22) ={22}; -Curve Loop(23) = { 29, 158, -30, -148}; -Plane Surface(23) ={23}; -Curve Loop(24) = { 30, -159, -31, 149}; -Plane Surface(24) ={24}; -Curve Loop(25) = { 31, 160, -32, -150}; -Plane Surface(25) ={25}; -Curve Loop(26) = { 32, -161, -33, 151}; -Plane Surface(26) ={26}; -Curve Loop(27) = { 33, 162, -34, -152}; -Plane Surface(27) ={27}; -Curve Loop(28) = { 34, -163, -35, 153}; -Plane Surface(28) ={28}; -Curve Loop(29) = { 35, 164, -36, -154}; -Plane Surface(29) ={29}; -Curve Loop(30) = { 36, -165, -37, 155}; -Plane Surface(30) ={30}; -Curve Loop(31) = { -38, 166, 39, -156}; -Plane Surface(31) ={31}; -Curve Loop(32) = { -39, -167, 40, 157}; -Plane Surface(32) ={32}; -Curve Loop(33) = { -40, 168, 41, -158}; -Plane Surface(33) ={33}; -Curve Loop(34) = { -41, -169, 42, 159}; -Plane Surface(34) ={34}; -Curve Loop(35) = { -42, 170, 43, -160}; -Plane Surface(35) ={35}; -Curve Loop(36) = { -43, -171, 44, 161}; -Plane Surface(36) ={36}; -Curve Loop(37) = { -44, 172, 45, -162}; -Plane Surface(37) ={37}; -Curve Loop(38) = { -45, -173, 46, 163}; -Plane Surface(38) ={38}; -Curve Loop(39) = { -46, 174, 47, -164}; -Plane Surface(39) ={39}; -Curve Loop(40) = { -47, -175, 48, 165}; -Plane Surface(40) ={40}; -Curve Loop(41) = { 49, 176, -50, -166}; -Plane Surface(41) ={41}; -Curve Loop(42) = { 50, -177, -51, 167}; -Plane Surface(42) ={42}; -Curve Loop(43) = { 51, 178, -52, -168}; -Plane Surface(43) ={43}; -Curve Loop(44) = { 52, -179, -53, 169}; -Plane Surface(44) ={44}; -Curve Loop(45) = { 53, 180, -54, -170}; -Plane Surface(45) ={45}; -Curve Loop(46) = { 54, -181, -55, 171}; -Plane Surface(46) ={46}; -Curve Loop(47) = { 55, 182, -56, -172}; -Plane Surface(47) ={47}; -Curve Loop(48) = { 56, -183, -57, 173}; -Plane Surface(48) ={48}; -Curve Loop(49) = { 57, 184, -58, -174}; -Plane Surface(49) ={49}; -Curve Loop(50) = { 58, -185, -59, 175}; -Plane Surface(50) ={50}; -Curve Loop(51) = { -60, 186, 61, -176}; -Plane Surface(51) ={51}; -Curve Loop(52) = { -61, -187, 62, 177}; -Plane Surface(52) ={52}; -Curve Loop(53) = { -62, 188, 63, -178}; -Plane Surface(53) ={53}; -Curve Loop(54) = { -63, -189, 64, 179}; -Plane Surface(54) ={54}; -Curve Loop(55) = { -64, 190, 65, -180}; -Plane Surface(55) ={55}; -Curve Loop(56) = { -65, -191, 66, 181}; -Plane Surface(56) ={56}; -Curve Loop(57) = { -66, 192, 67, -182}; -Plane Surface(57) ={57}; -Curve Loop(58) = { -67, -193, 68, 183}; -Plane Surface(58) ={58}; -Curve Loop(59) = { -68, 194, 69, -184}; -Plane Surface(59) ={59}; -Curve Loop(60) = { -69, -195, 70, 185}; -Plane Surface(60) ={60}; -Curve Loop(61) = { 71, 196, -72, -186}; -Plane Surface(61) ={61}; -Curve Loop(62) = { 72, -197, -73, 187}; -Plane Surface(62) ={62}; -Curve Loop(63) = { 73, 198, -74, -188}; -Plane Surface(63) ={63}; -Curve Loop(64) = { 74, -199, -75, 189}; -Plane Surface(64) ={64}; -Curve Loop(65) = { 75, 200, -76, -190}; -Plane Surface(65) ={65}; -Curve Loop(66) = { 76, -201, -77, 191}; -Plane Surface(66) ={66}; -Curve Loop(67) = { 77, 202, -78, -192}; -Plane Surface(67) ={67}; -Curve Loop(68) = { 78, -203, -79, 193}; -Plane Surface(68) ={68}; -Curve Loop(69) = { 79, 204, -80, -194}; -Plane Surface(69) ={69}; -Curve Loop(70) = { 80, -205, -81, 195}; -Plane Surface(70) ={70}; -Curve Loop(71) = { -82, 206, 83, -196}; -Plane Surface(71) ={71}; -Curve Loop(72) = { -83, -207, 84, 197}; -Plane Surface(72) ={72}; -Curve Loop(73) = { -84, 208, 85, -198}; -Plane Surface(73) ={73}; -Curve Loop(74) = { -85, -209, 86, 199}; -Plane Surface(74) ={74}; -Curve Loop(75) = { -86, 210, 87, -200}; -Plane Surface(75) ={75}; -Curve Loop(76) = { -87, -211, 88, 201}; -Plane Surface(76) ={76}; -Curve Loop(77) = { -88, 212, 89, -202}; -Plane Surface(77) ={77}; -Curve Loop(78) = { -89, -213, 90, 203}; -Plane Surface(78) ={78}; -Curve Loop(79) = { -90, 214, 91, -204}; -Plane Surface(79) ={79}; -Curve Loop(80) = { -91, -215, 92, 205}; -Plane Surface(80) ={80}; -Curve Loop(81) = { 93, 216, -94, -206}; -Plane Surface(81) ={81}; -Curve Loop(82) = { 94, -217, -95, 207}; -Plane Surface(82) ={82}; -Curve Loop(83) = { 95, 218, -96, -208}; -Plane Surface(83) ={83}; -Curve Loop(84) = { 96, -219, -97, 209}; -Plane Surface(84) ={84}; -Curve Loop(85) = { 97, 220, -98, -210}; -Plane Surface(85) ={85}; -Curve Loop(86) = { 98, -221, -99, 211}; -Plane Surface(86) ={86}; -Curve Loop(87) = { 99, 222, -100, -212}; -Plane Surface(87) ={87}; -Curve Loop(88) = { 100, -223, -101, 213}; -Plane Surface(88) ={88}; -Curve Loop(89) = { 101, 224, -102, -214}; -Plane Surface(89) ={89}; -Curve Loop(90) = { 102, -225, -103, 215}; -Plane Surface(90) ={90}; -Curve Loop(91) = { -104, 226, 105, -216}; -Plane Surface(91) ={91}; -Curve Loop(92) = { -105, -227, 106, 217}; -Plane Surface(92) ={92}; -Curve Loop(93) = { -106, 228, 107, -218}; -Plane Surface(93) ={93}; -Curve Loop(94) = { -107, -229, 108, 219}; -Plane Surface(94) ={94}; -Curve Loop(95) = { -108, 230, 109, -220}; -Plane Surface(95) ={95}; -Curve Loop(96) = { -109, -231, 110, 221}; -Plane Surface(96) ={96}; -Curve Loop(97) = { -110, 232, 111, -222}; -Plane Surface(97) ={97}; -Curve Loop(98) = { -111, -233, 112, 223}; -Plane Surface(98) ={98}; -Curve Loop(99) = { -112, 234, 113, -224}; -Plane Surface(99) ={99}; -Curve Loop(100) = { -113, -235, 114, 225}; -Plane Surface(100) ={100}; - - -// "Outside area" -//Curve Loop(101) = {1, 2, 3, 4}; -//+ -//Curve Loop(102) = {5, -16, 27, -38, 49, -60, 71, -82, 93, -104, 226, 227, 228, 229, 230, 231, 232, 233, 234, 235, 114, -103, 92, -81, 70, -59, 48, -37, 26, -15, -135, -134, -133, -132, -131, -130, -129, -128, -127, -126}; -//+ -//Plane Surface(101) = {101, 102}; - -//+ -Transfinite Surface {1}; -Transfinite Surface {2}; -Transfinite Surface {3}; -Transfinite Surface {4}; -Transfinite Surface {5}; -Transfinite Surface {6}; -Transfinite Surface {7}; -Transfinite Surface {8}; -Transfinite Surface {9}; -Transfinite Surface {10}; -Transfinite Surface {11}; -Transfinite Surface {12}; -Transfinite Surface {13}; -Transfinite Surface {14}; -Transfinite Surface {15}; -Transfinite Surface {16}; -Transfinite Surface {17}; -Transfinite Surface {18}; -Transfinite Surface {19}; -Transfinite Surface {20}; -Transfinite Surface {21}; -Transfinite Surface {22}; -Transfinite Surface {23}; -Transfinite Surface {24}; -Transfinite Surface {25}; -Transfinite Surface {26}; -Transfinite Surface {27}; -Transfinite Surface {28}; -Transfinite Surface {29}; -Transfinite Surface {30}; -Transfinite Surface {31}; -Transfinite Surface {32}; -Transfinite Surface {33}; -Transfinite Surface {34}; -Transfinite Surface {35}; -Transfinite Surface {36}; -Transfinite Surface {37}; -Transfinite Surface {38}; -Transfinite Surface {39}; -Transfinite Surface {40}; -Transfinite Surface {41}; -Transfinite Surface {42}; -Transfinite Surface {43}; -Transfinite Surface {44}; -Transfinite Surface {45}; -Transfinite Surface {46}; -Transfinite Surface {47}; -Transfinite Surface {48}; -Transfinite Surface {49}; -Transfinite Surface {50}; -Transfinite Surface {51}; -Transfinite Surface {52}; -Transfinite Surface {53}; -Transfinite Surface {54}; -Transfinite Surface {55}; -Transfinite Surface {56}; -Transfinite Surface {57}; -Transfinite Surface {58}; -Transfinite Surface {59}; -Transfinite Surface {60}; -Transfinite Surface {61}; -Transfinite Surface {62}; -Transfinite Surface {63}; -Transfinite Surface {64}; -Transfinite Surface {65}; -Transfinite Surface {66}; -Transfinite Surface {67}; -Transfinite Surface {68}; -Transfinite Surface {69}; -Transfinite Surface {70}; -Transfinite Surface {71}; -Transfinite Surface {72}; -Transfinite Surface {73}; -Transfinite Surface {74}; -Transfinite Surface {75}; -Transfinite Surface {76}; -Transfinite Surface {77}; -Transfinite Surface {78}; -Transfinite Surface {79}; -Transfinite Surface {80}; -Transfinite Surface {81}; -Transfinite Surface {82}; -Transfinite Surface {83}; -Transfinite Surface {84}; -Transfinite Surface {85}; -Transfinite Surface {86}; -Transfinite Surface {87}; -Transfinite Surface {88}; -Transfinite Surface {89}; -Transfinite Surface {90}; -Transfinite Surface {91}; -Transfinite Surface {92}; -Transfinite Surface {93}; -Transfinite Surface {94}; -Transfinite Surface {95}; -Transfinite Surface {96}; -Transfinite Surface {97}; -Transfinite Surface {98}; -Transfinite Surface {99}; -Transfinite Surface {100}; -//+ -Transfinite Curve {1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103} = n_recombine Using Progression n_prog; -// add outer perimeter with half as many cells and no Progression -Transfinite Curve {104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15} = n_recombine / 4 Using Progression 1; -//+ -Transfinite Curve { 127, 137, 147, 157, 167, 177, 187, 197, 207, 217, 227, 128, 138, 148, 158, 168, 178, 188, 198, 208, 218, 228, 129, 139, 149, 159, 169, 179, 189, 199, 209, 219, 229, 130, 140, 150, 160, 170, 180, 190, 200, 210, 220, 230, 131, 141, 151, 161, 171, 181, 191, 201, 211, 221, 231, 132, 142, 152, 162, 172, 182, 192, 202, 212, 222, 232, 133, 143, 153, 163, 173, 183, 193, 203, 213, 223, 233, 134, 144, 154, 164, 174, 184, 194, 204, 214, 224, 234} = n_recombine Using Progression n_prog; -// add outer perimeter with half as many cells and no Progression -Transfinite Curve {126, 136, 146, 156, 166, 176, 186, 196, 206, 216, 226, 135, 145, 155, 165, 175, 185, 195, 205, 215, 225, 235} = n_recombine / 4 Using Progression 1; -//+ -Line(236) = {1, 126}; -//+ -Line(237) = {126, 127}; -//+ -Line(238) = {127, 128}; -//+ -Line(239) = {128, 129}; -//+ -Line(240) = {129, 130}; -//+ -Line(241) = {130, 131}; -//+ -Line(242) = {131, 132}; -//+ -Line(243) = {132, 133}; -//+ -Line(244) = {133, 134}; -//+ -Line(245) = {134, 135}; -//+ -Line(246) = {135, 136}; -//+ -Line(247) = {136, 2}; -//+ -Line(248) = {2, 159}; -//+ -Line(249) = {159, 160}; -//+ -Line(250) = {160, 161}; -//+ -Line(251) = {161, 162}; -//+ -Line(252) = {162, 163}; -//+ -Line(253) = {163, 164}; -//+ -Line(254) = {164, 165}; -//+ -Line(255) = {165, 166}; -//+ -Line(256) = {166, 167}; -//+ -Line(257) = {167, 168}; -//+ -Line(258) = {168, 169}; -//+ -Line(259) = {169, 4}; -//+ -Line(260) = {4, 147}; -//+ -Line(261) = {147, 146}; -//+ -Line(262) = {146, 145}; -//+ -Line(263) = {145, 144}; -//+ -Line(264) = {144, 143}; -//+ -Line(265) = {143, 142}; -//+ -Line(266) = {142, 141}; -//+ -Line(267) = {141, 140}; -//+ -Line(268) = {140, 139}; -//+ -Line(269) = {139, 138}; -//+ -Line(270) = {138, 137}; -//+ -Line(271) = {137, 3}; -//+ -Line(272) = {3, 158}; -//+ -Line(273) = {158, 157}; -//+ -Line(274) = {157, 156}; -//+ -Line(275) = {156, 155}; -//+ -Line(276) = {155, 154}; -//+ -Line(277) = {154, 153}; -//+ -Line(278) = {153, 152}; -//+ -Line(279) = {152, 151}; -//+ -Line(280) = {151, 150}; -//+ -Line(281) = {150, 149}; -//+ -Line(282) = {149, 148}; -//+ -Line(283) = {148, 1}; -//+ -Line(284) = {126, 5}; -//+ -Line(285) = {127, 16}; -//+ -Line(286) = {128, 27}; -//+ -Line(287) = {129, 38}; -//+ -Line(288) = {130, 49}; -//+ -Line(289) = {131, 60}; -//+ -Line(290) = {132, 71}; -//+ -Line(291) = {133, 82}; -//+ -Line(292) = {134, 93}; -//+ -Line(293) = {135, 104}; -//+ -Line(294) = {136, 115}; -//+ -Line(295) = {159, 115}; -//+ -Line(296) = {160, 116}; -//+ -Line(297) = {161, 117}; -//+ -Line(298) = {162, 118}; -//+ -Line(299) = {163, 119}; -//+ -Line(300) = {164, 120}; -//+ -Line(301) = {165, 121}; -//+ -Line(302) = {166, 122}; -//+ -Line(303) = {167, 123}; -//+ -Line(304) = {168, 124}; -//+ -Line(305) = {169, 125}; -//+ -Line(306) = {147, 125}; -//+ -Line(307) = {146, 114}; -//+ -Line(308) = {145, 103}; -//+ -Line(309) = {144, 92}; -//+ -Line(310) = {143, 81}; -//+ -Line(311) = {142, 70}; -//+ -Line(312) = {141, 59}; -//+ -Line(313) = {140, 48}; -//+ -Line(314) = {139, 37}; -//+ -Line(315) = {138, 26}; -//+ -Line(316) = {137, 15}; -//+ -Line(317) = {158, 15}; -//+ -Line(318) = {157, 14}; -//+ -Line(319) = {156, 13}; -//+ -Line(320) = {155, 12}; -//+ -Line(321) = {154, 11}; -//+ -Line(322) = {153, 10}; -//+ -Line(323) = {152, 9}; -//+ -Line(324) = {151, 8}; -//+ -Line(325) = {150, 7}; -//+ -Line(326) = {149, 6}; -//+ -Line(327) = {148, 5}; -//+ -Curve Loop(103) = {236, 284, -327, 283}; -//+ -Plane Surface(102) = {103}; -//+ -Curve Loop(104) = {237, 285, -5, -284}; -//+ -Plane Surface(103) = {104}; -//+ -Curve Loop(105) = {238, 286, 16, -285}; -//+ -Plane Surface(104) = {105}; -//+ -Curve Loop(106) = {239, 287, -27, -286}; -//+ -Plane Surface(105) = {106}; -//+ -Curve Loop(107) = {240, 288, 38, -287}; -//+ -Plane Surface(106) = {107}; -//+ -Curve Loop(108) = {241, 289, -49, -288}; -//+ -Plane Surface(107) = {108}; -//+ -Curve Loop(109) = {242, 290, 60, -289}; -//+ -Plane Surface(108) = {109}; -//+ -Curve Loop(110) = {243, 291, -71, -290}; -//+ -Plane Surface(109) = {110}; -//+ -Curve Loop(111) = {244, 292, 82, -291}; -//+ -Plane Surface(110) = {111}; -//+ -Curve Loop(112) = {245, 293, -93, -292}; -//+ -Plane Surface(111) = {112}; -//+ -Curve Loop(113) = {246, 294, 104, -293}; -//+ -Plane Surface(112) = {113}; -//+ -Curve Loop(114) = {247, 248, 295, -294}; -//+ -Plane Surface(113) = {114}; -//+ -Curve Loop(115) = {249, 296, -226, -295}; -//+ -Plane Surface(114) = {115}; -//+ -Curve Loop(116) = {250, 297, 227, -296}; -//+ -Plane Surface(115) = {116}; -//+ -Curve Loop(117) = {251, 298, -228, -297}; -//+ -Plane Surface(116) = {117}; -//+ -Curve Loop(118) = {252, 299, 229, -298}; -//+ -Plane Surface(117) = {118}; -//+ -Curve Loop(119) = {253, 300, -230, -299}; -//+ -Plane Surface(118) = {119}; -//+ -Curve Loop(120) = {254, 301, 231, -300}; -//+ -Plane Surface(119) = {120}; -//+ -Curve Loop(121) = {255, 302, -232, -301}; -//+ -Plane Surface(120) = {121}; -//+ -Curve Loop(122) = {256, 303, 233, -302}; -//+ -Plane Surface(121) = {122}; -//+ -Curve Loop(123) = {257, 304, -234, -303}; -//+ -Plane Surface(122) = {123}; -//+ -Curve Loop(124) = {258, 305, 235, -304}; -//+ -Plane Surface(123) = {124}; -//+ -Curve Loop(125) = {259, 260, 306, -305}; -//+ -Plane Surface(124) = {125}; -//+ -Curve Loop(126) = {261, 307, -114, -306}; -//+ -Plane Surface(125) = {126}; -//+ -Curve Loop(127) = {262, 308, 103, -307}; -//+ -Plane Surface(126) = {127}; -//+ -Curve Loop(128) = {263, 309, -92, -308}; -//+ -Plane Surface(127) = {128}; -//+ -Curve Loop(129) = {264, 310, 81, -309}; -//+ -Plane Surface(128) = {129}; -//+ -Curve Loop(130) = {265, 311, -70, -310}; -//+ -Plane Surface(129) = {130}; -//+ -Curve Loop(131) = {266, 312, 59, -311}; -//+ -Plane Surface(130) = {131}; -//+ -Curve Loop(132) = {267, 313, -48, -312}; -//+ -Plane Surface(131) = {132}; -//+ -Curve Loop(133) = {268, 314, 37, -313}; -//+ -Plane Surface(132) = {133}; -//+ -Curve Loop(134) = {269, 315, -26, -314}; -//+ -Plane Surface(133) = {134}; -//+ -Curve Loop(135) = {270, 316, 15, -315}; -//+ -Plane Surface(134) = {135}; -//+ -Curve Loop(136) = {271, 272, 317, -316}; -//+ -Plane Surface(135) = {136}; -//+ -Curve Loop(137) = {273, 318, -135, -317}; -//+ -Plane Surface(136) = {137}; -//+ -Curve Loop(138) = {274, 319, 134, -318}; -//+ -Plane Surface(137) = {138}; -//+ -Curve Loop(139) = {275, 320, -133, -319}; -//+ -Plane Surface(138) = {139}; -//+ -Curve Loop(140) = {276, 321, 132, -320}; -//+ -Plane Surface(139) = {140}; -//+ -Curve Loop(141) = {277, 322, -131, -321}; -//+ -Plane Surface(140) = {141}; -//+ -Curve Loop(142) = {278, 323, 130, -322}; -//+ -Plane Surface(141) = {142}; -//+ -Curve Loop(143) = {279, 324, -129, -323}; -//+ -Plane Surface(142) = {143}; -//+ -Curve Loop(144) = {280, 325, 128, -324}; -//+ -Plane Surface(143) = {144}; -//+ -Curve Loop(145) = {281, 326, -127, -325}; -//+ -Plane Surface(144) = {145}; -//+ -Curve Loop(146) = {282, 327, 126, -326}; -//+ -Plane Surface(145) = {146}; -//+ -//+ -Transfinite Surface {102}; -//+ -Transfinite Surface {103}; -//+ -Transfinite Surface {104}; -//+ -Transfinite Surface {105}; -//+ -Transfinite Surface {106}; -//+ -Transfinite Surface {107}; -//+ -Transfinite Surface {108}; -//+ -Transfinite Surface {109}; -//+ -Transfinite Surface {110}; -//+ -Transfinite Surface {111}; -//+ -Transfinite Surface {112}; -//+ -Transfinite Surface {113}; -//+ -Transfinite Surface {114}; -//+ -Transfinite Surface {115}; -//+ -Transfinite Surface {116}; -//+ -Transfinite Surface {117}; -//+ -Transfinite Surface {118}; -//+ -Transfinite Surface {119}; -//+ -Transfinite Surface {120}; -//+ -Transfinite Surface {121}; -//+ -Transfinite Surface {122}; -//+ -Transfinite Surface {123}; -//+ -Transfinite Surface {124}; -//+ -Transfinite Surface {125}; -//+ -Transfinite Surface {126}; -//+ -Transfinite Surface {127}; -//+ -Transfinite Surface {128}; -//+ -Transfinite Surface {129}; -//+ -Transfinite Surface {130}; -//+ -Transfinite Surface {131}; -//+ -Transfinite Surface {132}; -//+ -Transfinite Surface {133}; -//+ -Transfinite Surface {134}; -//+ -Transfinite Surface {135}; -//+ -Transfinite Surface {136}; -//+ -Transfinite Surface {137}; -//+ -Transfinite Surface {138}; -//+ -Transfinite Surface {139}; -//+ -Transfinite Surface {140}; -//+ -Transfinite Surface {141}; -//+ -Transfinite Surface {142}; -//+ -Transfinite Surface {143}; -//+ -Transfinite Surface {144}; -//+ -Transfinite Surface {145}; -//+ -Transfinite Curve { 238, 239, 240, 241, 242, 243, 244, 245, 250, 251, 252, 253, 254, 255, 256, 257, 262, 263, 264, 265, 266, 267, 268, 269, 274, 275, 276, 277, 278, 279, 280, 281} = n_recombine Using Progression 1; -Transfinite Curve {237,246,249, 258,261,270,273,282} = n_recombine / 4 Using Progression 1; - - -Physical Curve("void",328) = {236, 237, 238, 239, 240, 241, 242, 243, 244, 245, 246, 247, 248, 249, 250, 251, 252, 253, 254, 255, 256, 257, 258, 259, 260, 261, 262, 263, 264, 265, 266, 267, 268, 269, 270, 271, 272, 273, 274, 275, 276, 277, 278, 279, 280, 281, 282, 283}; - -Recombine Surface "*"; -// Define meshing options -//Mesh.Algorithm = 6; // Specify meshing algorithm (e.g., Delaunay) -//Mesh.ElementOrder = 1; // Specify element order -// Generate the mesh -//Mesh 2; \ No newline at end of file diff --git a/examples_hpc/lattice/mesh/lattice_rectangular_n80.geo b/examples_hpc/lattice/mesh/lattice_rectangular_n80.geo deleted file mode 100644 index 19e33600..00000000 --- a/examples_hpc/lattice/mesh/lattice_rectangular_n80.geo +++ /dev/null @@ -1,1193 +0,0 @@ -cl_fine = 0.1; -n_recombine = 80; -n_prog = 1.05; -Point(1) = {-3.5, -3.5, 0, cl_fine}; -Point(2) = {3.5, -3.5, 0, cl_fine}; -Point(3) = {-3.5, 3.5, 0, cl_fine}; -Point(4) = {3.5, 3.5, 0, cl_fine}; - -// Inner grid -Point(5) = { -2.5, -2.5, 0, cl_fine }; -Point(6) = { -2.5, -2.0, 0, cl_fine }; -Point(7) = { -2.5, -1.5, 0, cl_fine }; -Point(8) = { -2.5, -1.0, 0, cl_fine }; -Point(9) = { -2.5, -0.5, 0, cl_fine }; -Point(10) = { -2.5, 0.0, 0, cl_fine }; -Point(11) = { -2.5, 0.5, 0, cl_fine }; -Point(12) = { -2.5, 1.0, 0, cl_fine }; -Point(13) = { -2.5, 1.5, 0, cl_fine }; -Point(14) = { -2.5, 2.0, 0, cl_fine }; -Point(15) = { -2.5, 2.5, 0, cl_fine }; -Point(16) = { -2.0, -2.5, 0, cl_fine }; -Point(17) = { -2.0, -2.0, 0, cl_fine }; -Point(18) = { -2.0, -1.5, 0, cl_fine }; -Point(19) = { -2.0, -1.0, 0, cl_fine }; -Point(20) = { -2.0, -0.5, 0, cl_fine }; -Point(21) = { -2.0, 0.0, 0, cl_fine }; -Point(22) = { -2.0, 0.5, 0, cl_fine }; -Point(23) = { -2.0, 1.0, 0, cl_fine }; -Point(24) = { -2.0, 1.5, 0, cl_fine }; -Point(25) = { -2.0, 2.0, 0, cl_fine }; -Point(26) = { -2.0, 2.5, 0, cl_fine }; -Point(27) = { -1.5, -2.5, 0, cl_fine }; -Point(28) = { -1.5, -2.0, 0, cl_fine }; -Point(29) = { -1.5, -1.5, 0, cl_fine }; -Point(30) = { -1.5, -1.0, 0, cl_fine }; -Point(31) = { -1.5, -0.5, 0, cl_fine }; -Point(32) = { -1.5, 0.0, 0, cl_fine }; -Point(33) = { -1.5, 0.5, 0, cl_fine }; -Point(34) = { -1.5, 1.0, 0, cl_fine }; -Point(35) = { -1.5, 1.5, 0, cl_fine }; -Point(36) = { -1.5, 2.0, 0, cl_fine }; -Point(37) = { -1.5, 2.5, 0, cl_fine }; -Point(38) = { -1.0, -2.5, 0, cl_fine }; -Point(39) = { -1.0, -2.0, 0, cl_fine }; -Point(40) = { -1.0, -1.5, 0, cl_fine }; -Point(41) = { -1.0, -1.0, 0, cl_fine }; -Point(42) = { -1.0, -0.5, 0, cl_fine }; -Point(43) = { -1.0, 0.0, 0, cl_fine }; -Point(44) = { -1.0, 0.5, 0, cl_fine }; -Point(45) = { -1.0, 1.0, 0, cl_fine }; -Point(46) = { -1.0, 1.5, 0, cl_fine }; -Point(47) = { -1.0, 2.0, 0, cl_fine }; -Point(48) = { -1.0, 2.5, 0, cl_fine }; -Point(49) = { -0.5, -2.5, 0, cl_fine }; -Point(50) = { -0.5, -2.0, 0, cl_fine }; -Point(51) = { -0.5, -1.5, 0, cl_fine }; -Point(52) = { -0.5, -1.0, 0, cl_fine }; -Point(53) = { -0.5, -0.5, 0, cl_fine }; -Point(54) = { -0.5, 0.0, 0, cl_fine }; -Point(55) = { -0.5, 0.5, 0, cl_fine }; -Point(56) = { -0.5, 1.0, 0, cl_fine }; -Point(57) = { -0.5, 1.5, 0, cl_fine }; -Point(58) = { -0.5, 2.0, 0, cl_fine }; -Point(59) = { -0.5, 2.5, 0, cl_fine }; -Point(60) = { 0.0, -2.5, 0, cl_fine }; -Point(61) = { 0.0, -2.0, 0, cl_fine }; -Point(62) = { 0.0, -1.5, 0, cl_fine }; -Point(63) = { 0.0, -1.0, 0, cl_fine }; -Point(64) = { 0.0, -0.5, 0, cl_fine }; -Point(65) = { 0.0, 0.0, 0, cl_fine }; -Point(66) = { 0.0, 0.5, 0, cl_fine }; -Point(67) = { 0.0, 1.0, 0, cl_fine }; -Point(68) = { 0.0, 1.5, 0, cl_fine }; -Point(69) = { 0.0, 2.0, 0, cl_fine }; -Point(70) = { 0.0, 2.5, 0, cl_fine }; -Point(71) = { 0.5, -2.5, 0, cl_fine }; -Point(72) = { 0.5, -2.0, 0, cl_fine }; -Point(73) = { 0.5, -1.5, 0, cl_fine }; -Point(74) = { 0.5, -1.0, 0, cl_fine }; -Point(75) = { 0.5, -0.5, 0, cl_fine }; -Point(76) = { 0.5, 0.0, 0, cl_fine }; -Point(77) = { 0.5, 0.5, 0, cl_fine }; -Point(78) = { 0.5, 1.0, 0, cl_fine }; -Point(79) = { 0.5, 1.5, 0, cl_fine }; -Point(80) = { 0.5, 2.0, 0, cl_fine }; -Point(81) = { 0.5, 2.5, 0, cl_fine }; -Point(82) = { 1.0, -2.5, 0, cl_fine }; -Point(83) = { 1.0, -2.0, 0, cl_fine }; -Point(84) = { 1.0, -1.5, 0, cl_fine }; -Point(85) = { 1.0, -1.0, 0, cl_fine }; -Point(86) = { 1.0, -0.5, 0, cl_fine }; -Point(87) = { 1.0, 0.0, 0, cl_fine }; -Point(88) = { 1.0, 0.5, 0, cl_fine }; -Point(89) = { 1.0, 1.0, 0, cl_fine }; -Point(90) = { 1.0, 1.5, 0, cl_fine }; -Point(91) = { 1.0, 2.0, 0, cl_fine }; -Point(92) = { 1.0, 2.5, 0, cl_fine }; -Point(93) = { 1.5, -2.5, 0, cl_fine }; -Point(94) = { 1.5, -2.0, 0, cl_fine }; -Point(95) = { 1.5, -1.5, 0, cl_fine }; -Point(96) = { 1.5, -1.0, 0, cl_fine }; -Point(97) = { 1.5, -0.5, 0, cl_fine }; -Point(98) = { 1.5, 0.0, 0, cl_fine }; -Point(99) = { 1.5, 0.5, 0, cl_fine }; -Point(100) = { 1.5, 1.0, 0, cl_fine }; -Point(101) = { 1.5, 1.5, 0, cl_fine }; -Point(102) = { 1.5, 2.0, 0, cl_fine }; -Point(103) = { 1.5, 2.5, 0, cl_fine }; -Point(104) = { 2.0, -2.5, 0, cl_fine }; -Point(105) = { 2.0, -2.0, 0, cl_fine }; -Point(106) = { 2.0, -1.5, 0, cl_fine }; -Point(107) = { 2.0, -1.0, 0, cl_fine }; -Point(108) = { 2.0, -0.5, 0, cl_fine }; -Point(109) = { 2.0, 0.0, 0, cl_fine }; -Point(110) = { 2.0, 0.5, 0, cl_fine }; -Point(111) = { 2.0, 1.0, 0, cl_fine }; -Point(112) = { 2.0, 1.5, 0, cl_fine }; -Point(113) = { 2.0, 2.0, 0, cl_fine }; -Point(114) = { 2.0, 2.5, 0, cl_fine }; -Point(115) = { 2.5, -2.5, 0, cl_fine }; -Point(116) = { 2.5, -2.0, 0, cl_fine }; -Point(117) = { 2.5, -1.5, 0, cl_fine }; -Point(118) = { 2.5, -1.0, 0, cl_fine }; -Point(119) = { 2.5, -0.5, 0, cl_fine }; -Point(120) = { 2.5, 0.0, 0, cl_fine }; -Point(121) = { 2.5, 0.5, 0, cl_fine }; -Point(122) = { 2.5, 1.0, 0, cl_fine }; -Point(123) = { 2.5, 1.5, 0, cl_fine }; -Point(124) = { 2.5, 2.0, 0, cl_fine }; -Point(125) = { 2.5, 2.5, 0, cl_fine }; - -// helper boundary points -Point(126) = { -2.5, -3.5, 0, cl_fine }; -Point(127) = { -2.0, -3.5, 0, cl_fine }; -Point(128) = { -1.5, -3.5, 0, cl_fine }; -Point(129) = { -1.0, -3.5, 0, cl_fine }; -Point(130) = { -0.5, -3.5, 0, cl_fine }; -Point(131) = { 0.0, -3.5, 0, cl_fine }; -Point(132) = { 0.5, -3.5, 0, cl_fine }; -Point(133) = { 1.0, -3.5, 0, cl_fine }; -Point(134) = { 1.5, -3.5, 0, cl_fine }; -Point(135) = { 2.0, -3.5, 0, cl_fine }; -Point(136) = { 2.5, -3.5, 0, cl_fine }; -Point(137) = { -2.5, 3.5, 0, cl_fine }; -Point(138) = { -2.0, 3.5, 0, cl_fine }; -Point(139) = { -1.5, 3.5, 0, cl_fine }; -Point(140) = { -1.0, 3.5, 0, cl_fine }; -Point(141) = { -0.5, 3.5, 0, cl_fine }; -Point(142) = { 0.0, 3.5, 0, cl_fine }; -Point(143) = { 0.5, 3.5, 0, cl_fine }; -Point(144) = { 1.0, 3.5, 0, cl_fine }; -Point(145) = { 1.5, 3.5, 0, cl_fine }; -Point(146) = { 2.0, 3.5, 0, cl_fine }; -Point(147) = { 2.5, 3.5, 0, cl_fine }; -Point(148) = { -3.5, -2.5, 0, cl_fine }; -Point(149) = { -3.5, -2.0, 0, cl_fine }; -Point(150) = { -3.5, -1.5, 0, cl_fine }; -Point(151) = { -3.5, -1.0, 0, cl_fine }; -Point(152) = { -3.5, -0.5, 0, cl_fine }; -Point(153) = { -3.5, 0.0, 0, cl_fine }; -Point(154) = { -3.5, 0.5, 0, cl_fine }; -Point(155) = { -3.5, 1.0, 0, cl_fine }; -Point(156) = { -3.5, 1.5, 0, cl_fine }; -Point(157) = { -3.5, 2.0, 0, cl_fine }; -Point(158) = { -3.5, 2.5, 0, cl_fine }; -Point(159) = { 3.5, -2.5, 0, cl_fine }; -Point(160) = { 3.5, -2.0, 0, cl_fine }; -Point(161) = { 3.5, -1.5, 0, cl_fine }; -Point(162) = { 3.5, -1.0, 0, cl_fine }; -Point(163) = { 3.5, -0.5, 0, cl_fine }; -Point(164) = { 3.5, 0.0, 0, cl_fine }; -Point(165) = { 3.5, 0.5, 0, cl_fine }; -Point(166) = { 3.5, 1.0, 0, cl_fine }; -Point(167) = { 3.5, 1.5, 0, cl_fine }; -Point(168) = { 3.5, 2.0, 0, cl_fine }; -Point(169) = { 3.5, 2.5, 0, cl_fine }; - -// Horizontal and vertical lines in inner grid -Line(5) = { 5, 16}; -Line(6) = { 6, 17}; -Line(7) = { 7, 18}; -Line(8) = { 8, 19}; -Line(9) = { 9, 20}; -Line(10) = { 10, 21}; -Line(11) = { 11, 22}; -Line(12) = { 12, 23}; -Line(13) = { 13, 24}; -Line(14) = { 14, 25}; -Line(15) = { 15, 26}; -Line(16) = { 27, 16}; -Line(17) = { 28, 17}; -Line(18) = { 29, 18}; -Line(19) = { 30, 19}; -Line(20) = { 31, 20}; -Line(21) = { 32, 21}; -Line(22) = { 33, 22}; -Line(23) = { 34, 23}; -Line(24) = { 35, 24}; -Line(25) = { 36, 25}; -Line(26) = { 37, 26}; -Line(27) = { 27, 38}; -Line(28) = { 28, 39}; -Line(29) = { 29, 40}; -Line(30) = { 30, 41}; -Line(31) = { 31, 42}; -Line(32) = { 32, 43}; -Line(33) = { 33, 44}; -Line(34) = { 34, 45}; -Line(35) = { 35, 46}; -Line(36) = { 36, 47}; -Line(37) = { 37, 48}; -Line(38) = { 49, 38}; -Line(39) = { 50, 39}; -Line(40) = { 51, 40}; -Line(41) = { 52, 41}; -Line(42) = { 53, 42}; -Line(43) = { 54, 43}; -Line(44) = { 55, 44}; -Line(45) = { 56, 45}; -Line(46) = { 57, 46}; -Line(47) = { 58, 47}; -Line(48) = { 59, 48}; -Line(49) = { 49, 60}; -Line(50) = { 50, 61}; -Line(51) = { 51, 62}; -Line(52) = { 52, 63}; -Line(53) = { 53, 64}; -Line(54) = { 54, 65}; -Line(55) = { 55, 66}; -Line(56) = { 56, 67}; -Line(57) = { 57, 68}; -Line(58) = { 58, 69}; -Line(59) = { 59, 70}; -Line(60) = { 71, 60}; -Line(61) = { 72, 61}; -Line(62) = { 73, 62}; -Line(63) = { 74, 63}; -Line(64) = { 75, 64}; -Line(65) = { 76, 65}; -Line(66) = { 77, 66}; -Line(67) = { 78, 67}; -Line(68) = { 79, 68}; -Line(69) = { 80, 69}; -Line(70) = { 81, 70}; -Line(71) = { 71, 82}; -Line(72) = { 72, 83}; -Line(73) = { 73, 84}; -Line(74) = { 74, 85}; -Line(75) = { 75, 86}; -Line(76) = { 76, 87}; -Line(77) = { 77, 88}; -Line(78) = { 78, 89}; -Line(79) = { 79, 90}; -Line(80) = { 80, 91}; -Line(81) = { 81, 92}; -Line(82) = { 93, 82}; -Line(83) = { 94, 83}; -Line(84) = { 95, 84}; -Line(85) = { 96, 85}; -Line(86) = { 97, 86}; -Line(87) = { 98, 87}; -Line(88) = { 99, 88}; -Line(89) = { 100, 89}; -Line(90) = { 101, 90}; -Line(91) = { 102, 91}; -Line(92) = { 103, 92}; -Line(93) = { 93, 104}; -Line(94) = { 94, 105}; -Line(95) = { 95, 106}; -Line(96) = { 96, 107}; -Line(97) = { 97, 108}; -Line(98) = { 98, 109}; -Line(99) = { 99, 110}; -Line(100) = { 100, 111}; -Line(101) = { 101, 112}; -Line(102) = { 102, 113}; -Line(103) = { 103, 114}; -Line(104) = { 115, 104}; -Line(105) = { 116, 105}; -Line(106) = { 117, 106}; -Line(107) = { 118, 107}; -Line(108) = { 119, 108}; -Line(109) = { 120, 109}; -Line(110) = { 121, 110}; -Line(111) = { 122, 111}; -Line(112) = { 123, 112}; -Line(113) = { 124, 113}; -Line(114) = { 125, 114}; -Line(126) = { 5, 6}; -Line(127) = { 7, 6}; -Line(128) = { 7, 8}; -Line(129) = { 9, 8}; -Line(130) = { 9, 10}; -Line(131) = { 11, 10}; -Line(132) = { 11, 12}; -Line(133) = { 13, 12}; -Line(134) = { 13, 14}; -Line(135) = { 15, 14}; -Line(136) = { 16, 17}; -Line(137) = { 18, 17}; -Line(138) = { 18, 19}; -Line(139) = { 20, 19}; -Line(140) = { 20, 21}; -Line(141) = { 22, 21}; -Line(142) = { 22, 23}; -Line(143) = { 24, 23}; -Line(144) = { 24, 25}; -Line(145) = { 26, 25}; - -Line(146) = { 27, 28}; -Line(147) = { 29, 28}; -Line(148) = { 29, 30}; -Line(149) = { 31, 30}; -Line(150) = { 31, 32}; -Line(151) = { 33, 32}; -Line(152) = { 33, 34}; -Line(153) = { 35, 34}; -Line(154) = { 35, 36}; -Line(155) = { 37, 36}; - -Line(156) = { 38, 39}; -Line(157) = { 40, 39}; -Line(158) = { 40, 41}; -Line(159) = { 42, 41}; -Line(160) = { 42, 43}; -Line(161) = { 44, 43}; -Line(162) = { 44, 45}; -Line(163) = { 46, 45}; -Line(164) = { 46, 47}; -Line(165) = { 48, 47}; - -Line(166) = { 49, 50}; -Line(167) = { 51, 50}; -Line(168) = { 51, 52}; -Line(169) = { 53, 52}; -Line(170) = { 53, 54}; -Line(171) = { 55, 54}; -Line(172) = { 55, 56}; -Line(173) = { 57, 56}; -Line(174) = { 57, 58}; -Line(175) = { 59, 58}; - -Line(176) = { 60, 61}; -Line(177) = { 62, 61}; -Line(178) = { 62, 63}; -Line(179) = { 64, 63}; -Line(180) = { 64, 65}; -Line(181) = { 66, 65}; -Line(182) = { 66, 67}; -Line(183) = { 68, 67}; -Line(184) = { 68, 69}; -Line(185) = { 70, 69}; - -Line(186) = { 71, 72}; -Line(187) = { 73, 72}; -Line(188) = { 73, 74}; -Line(189) = { 75, 74}; -Line(190) = { 75, 76}; -Line(191) = { 77, 76}; -Line(192) = { 77, 78}; -Line(193) = { 79, 78}; -Line(194) = { 79, 80}; -Line(195) = { 81, 80}; - -Line(196) = { 82, 83}; -Line(197) = { 84, 83}; -Line(198) = { 84, 85}; -Line(199) = { 86, 85}; -Line(200) = { 86, 87}; -Line(201) = { 88, 87}; -Line(202) = { 88, 89}; -Line(203) = { 90, 89}; -Line(204) = { 90, 91}; -Line(205) = { 92, 91}; - -Line(206) = { 93, 94}; -Line(207) = { 95, 94}; -Line(208) = { 95, 96}; -Line(209) = { 97, 96}; -Line(210) = { 97, 98}; -Line(211) = { 99, 98}; -Line(212) = { 99, 100}; -Line(213) = { 101, 100}; -Line(214) = { 101, 102}; -Line(215) = { 103, 102}; - - -Line(216) = { 104, 105}; -Line(217) = { 106, 105}; -Line(218) = { 106, 107}; -Line(219) = { 108, 107}; -Line(220) = { 108, 109}; -Line(221) = { 110, 109}; -Line(222) = { 110, 111}; -Line(223) = { 112, 111}; -Line(224) = { 112, 113}; -Line(225) = { 114, 113}; - - -Line(226) = { 115, 116}; -Line(227) = { 117, 116}; -Line(228) = { 117, 118}; -Line(229) = { 119, 118}; -Line(230) = { 119, 120}; -Line(231) = { 121, 120}; -Line(232) = { 121, 122}; -Line(233) = { 123, 122}; -Line(234) = { 123, 124}; -Line(235) = { 125, 124}; - - -// Define inner surfaces -//+ -Curve Loop(1) = { 5, 136, -6, -126}; -Plane Surface(1) ={1}; -Curve Loop(2) = { 6, -137, -7, 127}; -Plane Surface(2) ={2}; -Curve Loop(3) = { 7, 138, -8, -128}; -Plane Surface(3) ={3}; -Curve Loop(4) = { 8, -139, -9, 129}; -Plane Surface(4) ={4}; -Curve Loop(5) = { 9, 140, -10, -130}; -Plane Surface(5) ={5}; -Curve Loop(6) = { 10, -141, -11, 131}; -Plane Surface(6) ={6}; -Curve Loop(7) = { 11, 142, -12, -132}; -Plane Surface(7) ={7}; -Curve Loop(8) = { 12, -143, -13, 133}; -Plane Surface(8) ={8}; -Curve Loop(9) = { 13, 144, -14, -134}; -Plane Surface(9) ={9}; -Curve Loop(10) = { 14, -145, -15, 135}; -Plane Surface(10) ={10}; -Curve Loop(11) = { -16, 146, 17, -136}; -Plane Surface(11) ={11}; -Curve Loop(12) = { -17, -147, 18, 137}; -Plane Surface(12) ={12}; -Curve Loop(13) = { -18, 148, 19, -138}; -Plane Surface(13) ={13}; -Curve Loop(14) = { -19, -149, 20, 139}; -Plane Surface(14) ={14}; -Curve Loop(15) = { -20, 150, 21, -140}; -Plane Surface(15) ={15}; -Curve Loop(16) = { -21, -151, 22, 141}; -Plane Surface(16) ={16}; -Curve Loop(17) = { -22, 152, 23, -142}; -Plane Surface(17) ={17}; -Curve Loop(18) = { -23, -153, 24, 143}; -Plane Surface(18) ={18}; -Curve Loop(19) = {-24, 154, 25, -144}; -Plane Surface(19) ={19}; -Curve Loop(20) = { -25, -155, 26, 145}; -Plane Surface(20) ={20}; -Curve Loop(21) = { 27, 156, -28, -146}; -Plane Surface(21) ={21}; -Curve Loop(22) = { 28, -157, -29, 147}; -Plane Surface(22) ={22}; -Curve Loop(23) = { 29, 158, -30, -148}; -Plane Surface(23) ={23}; -Curve Loop(24) = { 30, -159, -31, 149}; -Plane Surface(24) ={24}; -Curve Loop(25) = { 31, 160, -32, -150}; -Plane Surface(25) ={25}; -Curve Loop(26) = { 32, -161, -33, 151}; -Plane Surface(26) ={26}; -Curve Loop(27) = { 33, 162, -34, -152}; -Plane Surface(27) ={27}; -Curve Loop(28) = { 34, -163, -35, 153}; -Plane Surface(28) ={28}; -Curve Loop(29) = { 35, 164, -36, -154}; -Plane Surface(29) ={29}; -Curve Loop(30) = { 36, -165, -37, 155}; -Plane Surface(30) ={30}; -Curve Loop(31) = { -38, 166, 39, -156}; -Plane Surface(31) ={31}; -Curve Loop(32) = { -39, -167, 40, 157}; -Plane Surface(32) ={32}; -Curve Loop(33) = { -40, 168, 41, -158}; -Plane Surface(33) ={33}; -Curve Loop(34) = { -41, -169, 42, 159}; -Plane Surface(34) ={34}; -Curve Loop(35) = { -42, 170, 43, -160}; -Plane Surface(35) ={35}; -Curve Loop(36) = { -43, -171, 44, 161}; -Plane Surface(36) ={36}; -Curve Loop(37) = { -44, 172, 45, -162}; -Plane Surface(37) ={37}; -Curve Loop(38) = { -45, -173, 46, 163}; -Plane Surface(38) ={38}; -Curve Loop(39) = { -46, 174, 47, -164}; -Plane Surface(39) ={39}; -Curve Loop(40) = { -47, -175, 48, 165}; -Plane Surface(40) ={40}; -Curve Loop(41) = { 49, 176, -50, -166}; -Plane Surface(41) ={41}; -Curve Loop(42) = { 50, -177, -51, 167}; -Plane Surface(42) ={42}; -Curve Loop(43) = { 51, 178, -52, -168}; -Plane Surface(43) ={43}; -Curve Loop(44) = { 52, -179, -53, 169}; -Plane Surface(44) ={44}; -Curve Loop(45) = { 53, 180, -54, -170}; -Plane Surface(45) ={45}; -Curve Loop(46) = { 54, -181, -55, 171}; -Plane Surface(46) ={46}; -Curve Loop(47) = { 55, 182, -56, -172}; -Plane Surface(47) ={47}; -Curve Loop(48) = { 56, -183, -57, 173}; -Plane Surface(48) ={48}; -Curve Loop(49) = { 57, 184, -58, -174}; -Plane Surface(49) ={49}; -Curve Loop(50) = { 58, -185, -59, 175}; -Plane Surface(50) ={50}; -Curve Loop(51) = { -60, 186, 61, -176}; -Plane Surface(51) ={51}; -Curve Loop(52) = { -61, -187, 62, 177}; -Plane Surface(52) ={52}; -Curve Loop(53) = { -62, 188, 63, -178}; -Plane Surface(53) ={53}; -Curve Loop(54) = { -63, -189, 64, 179}; -Plane Surface(54) ={54}; -Curve Loop(55) = { -64, 190, 65, -180}; -Plane Surface(55) ={55}; -Curve Loop(56) = { -65, -191, 66, 181}; -Plane Surface(56) ={56}; -Curve Loop(57) = { -66, 192, 67, -182}; -Plane Surface(57) ={57}; -Curve Loop(58) = { -67, -193, 68, 183}; -Plane Surface(58) ={58}; -Curve Loop(59) = { -68, 194, 69, -184}; -Plane Surface(59) ={59}; -Curve Loop(60) = { -69, -195, 70, 185}; -Plane Surface(60) ={60}; -Curve Loop(61) = { 71, 196, -72, -186}; -Plane Surface(61) ={61}; -Curve Loop(62) = { 72, -197, -73, 187}; -Plane Surface(62) ={62}; -Curve Loop(63) = { 73, 198, -74, -188}; -Plane Surface(63) ={63}; -Curve Loop(64) = { 74, -199, -75, 189}; -Plane Surface(64) ={64}; -Curve Loop(65) = { 75, 200, -76, -190}; -Plane Surface(65) ={65}; -Curve Loop(66) = { 76, -201, -77, 191}; -Plane Surface(66) ={66}; -Curve Loop(67) = { 77, 202, -78, -192}; -Plane Surface(67) ={67}; -Curve Loop(68) = { 78, -203, -79, 193}; -Plane Surface(68) ={68}; -Curve Loop(69) = { 79, 204, -80, -194}; -Plane Surface(69) ={69}; -Curve Loop(70) = { 80, -205, -81, 195}; -Plane Surface(70) ={70}; -Curve Loop(71) = { -82, 206, 83, -196}; -Plane Surface(71) ={71}; -Curve Loop(72) = { -83, -207, 84, 197}; -Plane Surface(72) ={72}; -Curve Loop(73) = { -84, 208, 85, -198}; -Plane Surface(73) ={73}; -Curve Loop(74) = { -85, -209, 86, 199}; -Plane Surface(74) ={74}; -Curve Loop(75) = { -86, 210, 87, -200}; -Plane Surface(75) ={75}; -Curve Loop(76) = { -87, -211, 88, 201}; -Plane Surface(76) ={76}; -Curve Loop(77) = { -88, 212, 89, -202}; -Plane Surface(77) ={77}; -Curve Loop(78) = { -89, -213, 90, 203}; -Plane Surface(78) ={78}; -Curve Loop(79) = { -90, 214, 91, -204}; -Plane Surface(79) ={79}; -Curve Loop(80) = { -91, -215, 92, 205}; -Plane Surface(80) ={80}; -Curve Loop(81) = { 93, 216, -94, -206}; -Plane Surface(81) ={81}; -Curve Loop(82) = { 94, -217, -95, 207}; -Plane Surface(82) ={82}; -Curve Loop(83) = { 95, 218, -96, -208}; -Plane Surface(83) ={83}; -Curve Loop(84) = { 96, -219, -97, 209}; -Plane Surface(84) ={84}; -Curve Loop(85) = { 97, 220, -98, -210}; -Plane Surface(85) ={85}; -Curve Loop(86) = { 98, -221, -99, 211}; -Plane Surface(86) ={86}; -Curve Loop(87) = { 99, 222, -100, -212}; -Plane Surface(87) ={87}; -Curve Loop(88) = { 100, -223, -101, 213}; -Plane Surface(88) ={88}; -Curve Loop(89) = { 101, 224, -102, -214}; -Plane Surface(89) ={89}; -Curve Loop(90) = { 102, -225, -103, 215}; -Plane Surface(90) ={90}; -Curve Loop(91) = { -104, 226, 105, -216}; -Plane Surface(91) ={91}; -Curve Loop(92) = { -105, -227, 106, 217}; -Plane Surface(92) ={92}; -Curve Loop(93) = { -106, 228, 107, -218}; -Plane Surface(93) ={93}; -Curve Loop(94) = { -107, -229, 108, 219}; -Plane Surface(94) ={94}; -Curve Loop(95) = { -108, 230, 109, -220}; -Plane Surface(95) ={95}; -Curve Loop(96) = { -109, -231, 110, 221}; -Plane Surface(96) ={96}; -Curve Loop(97) = { -110, 232, 111, -222}; -Plane Surface(97) ={97}; -Curve Loop(98) = { -111, -233, 112, 223}; -Plane Surface(98) ={98}; -Curve Loop(99) = { -112, 234, 113, -224}; -Plane Surface(99) ={99}; -Curve Loop(100) = { -113, -235, 114, 225}; -Plane Surface(100) ={100}; - - -// "Outside area" -//Curve Loop(101) = {1, 2, 3, 4}; -//+ -//Curve Loop(102) = {5, -16, 27, -38, 49, -60, 71, -82, 93, -104, 226, 227, 228, 229, 230, 231, 232, 233, 234, 235, 114, -103, 92, -81, 70, -59, 48, -37, 26, -15, -135, -134, -133, -132, -131, -130, -129, -128, -127, -126}; -//+ -//Plane Surface(101) = {101, 102}; - -//+ -Transfinite Surface {1}; -Transfinite Surface {2}; -Transfinite Surface {3}; -Transfinite Surface {4}; -Transfinite Surface {5}; -Transfinite Surface {6}; -Transfinite Surface {7}; -Transfinite Surface {8}; -Transfinite Surface {9}; -Transfinite Surface {10}; -Transfinite Surface {11}; -Transfinite Surface {12}; -Transfinite Surface {13}; -Transfinite Surface {14}; -Transfinite Surface {15}; -Transfinite Surface {16}; -Transfinite Surface {17}; -Transfinite Surface {18}; -Transfinite Surface {19}; -Transfinite Surface {20}; -Transfinite Surface {21}; -Transfinite Surface {22}; -Transfinite Surface {23}; -Transfinite Surface {24}; -Transfinite Surface {25}; -Transfinite Surface {26}; -Transfinite Surface {27}; -Transfinite Surface {28}; -Transfinite Surface {29}; -Transfinite Surface {30}; -Transfinite Surface {31}; -Transfinite Surface {32}; -Transfinite Surface {33}; -Transfinite Surface {34}; -Transfinite Surface {35}; -Transfinite Surface {36}; -Transfinite Surface {37}; -Transfinite Surface {38}; -Transfinite Surface {39}; -Transfinite Surface {40}; -Transfinite Surface {41}; -Transfinite Surface {42}; -Transfinite Surface {43}; -Transfinite Surface {44}; -Transfinite Surface {45}; -Transfinite Surface {46}; -Transfinite Surface {47}; -Transfinite Surface {48}; -Transfinite Surface {49}; -Transfinite Surface {50}; -Transfinite Surface {51}; -Transfinite Surface {52}; -Transfinite Surface {53}; -Transfinite Surface {54}; -Transfinite Surface {55}; -Transfinite Surface {56}; -Transfinite Surface {57}; -Transfinite Surface {58}; -Transfinite Surface {59}; -Transfinite Surface {60}; -Transfinite Surface {61}; -Transfinite Surface {62}; -Transfinite Surface {63}; -Transfinite Surface {64}; -Transfinite Surface {65}; -Transfinite Surface {66}; -Transfinite Surface {67}; -Transfinite Surface {68}; -Transfinite Surface {69}; -Transfinite Surface {70}; -Transfinite Surface {71}; -Transfinite Surface {72}; -Transfinite Surface {73}; -Transfinite Surface {74}; -Transfinite Surface {75}; -Transfinite Surface {76}; -Transfinite Surface {77}; -Transfinite Surface {78}; -Transfinite Surface {79}; -Transfinite Surface {80}; -Transfinite Surface {81}; -Transfinite Surface {82}; -Transfinite Surface {83}; -Transfinite Surface {84}; -Transfinite Surface {85}; -Transfinite Surface {86}; -Transfinite Surface {87}; -Transfinite Surface {88}; -Transfinite Surface {89}; -Transfinite Surface {90}; -Transfinite Surface {91}; -Transfinite Surface {92}; -Transfinite Surface {93}; -Transfinite Surface {94}; -Transfinite Surface {95}; -Transfinite Surface {96}; -Transfinite Surface {97}; -Transfinite Surface {98}; -Transfinite Surface {99}; -Transfinite Surface {100}; -//+ -Transfinite Curve {1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103} = n_recombine Using Progression n_prog; -// add outer perimeter with half as many cells and no Progression -Transfinite Curve {104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15} = n_recombine / 4 Using Progression 1; -//+ -Transfinite Curve { 127, 137, 147, 157, 167, 177, 187, 197, 207, 217, 227, 128, 138, 148, 158, 168, 178, 188, 198, 208, 218, 228, 129, 139, 149, 159, 169, 179, 189, 199, 209, 219, 229, 130, 140, 150, 160, 170, 180, 190, 200, 210, 220, 230, 131, 141, 151, 161, 171, 181, 191, 201, 211, 221, 231, 132, 142, 152, 162, 172, 182, 192, 202, 212, 222, 232, 133, 143, 153, 163, 173, 183, 193, 203, 213, 223, 233, 134, 144, 154, 164, 174, 184, 194, 204, 214, 224, 234} = n_recombine Using Progression n_prog; -// add outer perimeter with half as many cells and no Progression -Transfinite Curve {126, 136, 146, 156, 166, 176, 186, 196, 206, 216, 226, 135, 145, 155, 165, 175, 185, 195, 205, 215, 225, 235} = n_recombine / 4 Using Progression 1; -//+ -Line(236) = {1, 126}; -//+ -Line(237) = {126, 127}; -//+ -Line(238) = {127, 128}; -//+ -Line(239) = {128, 129}; -//+ -Line(240) = {129, 130}; -//+ -Line(241) = {130, 131}; -//+ -Line(242) = {131, 132}; -//+ -Line(243) = {132, 133}; -//+ -Line(244) = {133, 134}; -//+ -Line(245) = {134, 135}; -//+ -Line(246) = {135, 136}; -//+ -Line(247) = {136, 2}; -//+ -Line(248) = {2, 159}; -//+ -Line(249) = {159, 160}; -//+ -Line(250) = {160, 161}; -//+ -Line(251) = {161, 162}; -//+ -Line(252) = {162, 163}; -//+ -Line(253) = {163, 164}; -//+ -Line(254) = {164, 165}; -//+ -Line(255) = {165, 166}; -//+ -Line(256) = {166, 167}; -//+ -Line(257) = {167, 168}; -//+ -Line(258) = {168, 169}; -//+ -Line(259) = {169, 4}; -//+ -Line(260) = {4, 147}; -//+ -Line(261) = {147, 146}; -//+ -Line(262) = {146, 145}; -//+ -Line(263) = {145, 144}; -//+ -Line(264) = {144, 143}; -//+ -Line(265) = {143, 142}; -//+ -Line(266) = {142, 141}; -//+ -Line(267) = {141, 140}; -//+ -Line(268) = {140, 139}; -//+ -Line(269) = {139, 138}; -//+ -Line(270) = {138, 137}; -//+ -Line(271) = {137, 3}; -//+ -Line(272) = {3, 158}; -//+ -Line(273) = {158, 157}; -//+ -Line(274) = {157, 156}; -//+ -Line(275) = {156, 155}; -//+ -Line(276) = {155, 154}; -//+ -Line(277) = {154, 153}; -//+ -Line(278) = {153, 152}; -//+ -Line(279) = {152, 151}; -//+ -Line(280) = {151, 150}; -//+ -Line(281) = {150, 149}; -//+ -Line(282) = {149, 148}; -//+ -Line(283) = {148, 1}; -//+ -Line(284) = {126, 5}; -//+ -Line(285) = {127, 16}; -//+ -Line(286) = {128, 27}; -//+ -Line(287) = {129, 38}; -//+ -Line(288) = {130, 49}; -//+ -Line(289) = {131, 60}; -//+ -Line(290) = {132, 71}; -//+ -Line(291) = {133, 82}; -//+ -Line(292) = {134, 93}; -//+ -Line(293) = {135, 104}; -//+ -Line(294) = {136, 115}; -//+ -Line(295) = {159, 115}; -//+ -Line(296) = {160, 116}; -//+ -Line(297) = {161, 117}; -//+ -Line(298) = {162, 118}; -//+ -Line(299) = {163, 119}; -//+ -Line(300) = {164, 120}; -//+ -Line(301) = {165, 121}; -//+ -Line(302) = {166, 122}; -//+ -Line(303) = {167, 123}; -//+ -Line(304) = {168, 124}; -//+ -Line(305) = {169, 125}; -//+ -Line(306) = {147, 125}; -//+ -Line(307) = {146, 114}; -//+ -Line(308) = {145, 103}; -//+ -Line(309) = {144, 92}; -//+ -Line(310) = {143, 81}; -//+ -Line(311) = {142, 70}; -//+ -Line(312) = {141, 59}; -//+ -Line(313) = {140, 48}; -//+ -Line(314) = {139, 37}; -//+ -Line(315) = {138, 26}; -//+ -Line(316) = {137, 15}; -//+ -Line(317) = {158, 15}; -//+ -Line(318) = {157, 14}; -//+ -Line(319) = {156, 13}; -//+ -Line(320) = {155, 12}; -//+ -Line(321) = {154, 11}; -//+ -Line(322) = {153, 10}; -//+ -Line(323) = {152, 9}; -//+ -Line(324) = {151, 8}; -//+ -Line(325) = {150, 7}; -//+ -Line(326) = {149, 6}; -//+ -Line(327) = {148, 5}; -//+ -Curve Loop(103) = {236, 284, -327, 283}; -//+ -Plane Surface(102) = {103}; -//+ -Curve Loop(104) = {237, 285, -5, -284}; -//+ -Plane Surface(103) = {104}; -//+ -Curve Loop(105) = {238, 286, 16, -285}; -//+ -Plane Surface(104) = {105}; -//+ -Curve Loop(106) = {239, 287, -27, -286}; -//+ -Plane Surface(105) = {106}; -//+ -Curve Loop(107) = {240, 288, 38, -287}; -//+ -Plane Surface(106) = {107}; -//+ -Curve Loop(108) = {241, 289, -49, -288}; -//+ -Plane Surface(107) = {108}; -//+ -Curve Loop(109) = {242, 290, 60, -289}; -//+ -Plane Surface(108) = {109}; -//+ -Curve Loop(110) = {243, 291, -71, -290}; -//+ -Plane Surface(109) = {110}; -//+ -Curve Loop(111) = {244, 292, 82, -291}; -//+ -Plane Surface(110) = {111}; -//+ -Curve Loop(112) = {245, 293, -93, -292}; -//+ -Plane Surface(111) = {112}; -//+ -Curve Loop(113) = {246, 294, 104, -293}; -//+ -Plane Surface(112) = {113}; -//+ -Curve Loop(114) = {247, 248, 295, -294}; -//+ -Plane Surface(113) = {114}; -//+ -Curve Loop(115) = {249, 296, -226, -295}; -//+ -Plane Surface(114) = {115}; -//+ -Curve Loop(116) = {250, 297, 227, -296}; -//+ -Plane Surface(115) = {116}; -//+ -Curve Loop(117) = {251, 298, -228, -297}; -//+ -Plane Surface(116) = {117}; -//+ -Curve Loop(118) = {252, 299, 229, -298}; -//+ -Plane Surface(117) = {118}; -//+ -Curve Loop(119) = {253, 300, -230, -299}; -//+ -Plane Surface(118) = {119}; -//+ -Curve Loop(120) = {254, 301, 231, -300}; -//+ -Plane Surface(119) = {120}; -//+ -Curve Loop(121) = {255, 302, -232, -301}; -//+ -Plane Surface(120) = {121}; -//+ -Curve Loop(122) = {256, 303, 233, -302}; -//+ -Plane Surface(121) = {122}; -//+ -Curve Loop(123) = {257, 304, -234, -303}; -//+ -Plane Surface(122) = {123}; -//+ -Curve Loop(124) = {258, 305, 235, -304}; -//+ -Plane Surface(123) = {124}; -//+ -Curve Loop(125) = {259, 260, 306, -305}; -//+ -Plane Surface(124) = {125}; -//+ -Curve Loop(126) = {261, 307, -114, -306}; -//+ -Plane Surface(125) = {126}; -//+ -Curve Loop(127) = {262, 308, 103, -307}; -//+ -Plane Surface(126) = {127}; -//+ -Curve Loop(128) = {263, 309, -92, -308}; -//+ -Plane Surface(127) = {128}; -//+ -Curve Loop(129) = {264, 310, 81, -309}; -//+ -Plane Surface(128) = {129}; -//+ -Curve Loop(130) = {265, 311, -70, -310}; -//+ -Plane Surface(129) = {130}; -//+ -Curve Loop(131) = {266, 312, 59, -311}; -//+ -Plane Surface(130) = {131}; -//+ -Curve Loop(132) = {267, 313, -48, -312}; -//+ -Plane Surface(131) = {132}; -//+ -Curve Loop(133) = {268, 314, 37, -313}; -//+ -Plane Surface(132) = {133}; -//+ -Curve Loop(134) = {269, 315, -26, -314}; -//+ -Plane Surface(133) = {134}; -//+ -Curve Loop(135) = {270, 316, 15, -315}; -//+ -Plane Surface(134) = {135}; -//+ -Curve Loop(136) = {271, 272, 317, -316}; -//+ -Plane Surface(135) = {136}; -//+ -Curve Loop(137) = {273, 318, -135, -317}; -//+ -Plane Surface(136) = {137}; -//+ -Curve Loop(138) = {274, 319, 134, -318}; -//+ -Plane Surface(137) = {138}; -//+ -Curve Loop(139) = {275, 320, -133, -319}; -//+ -Plane Surface(138) = {139}; -//+ -Curve Loop(140) = {276, 321, 132, -320}; -//+ -Plane Surface(139) = {140}; -//+ -Curve Loop(141) = {277, 322, -131, -321}; -//+ -Plane Surface(140) = {141}; -//+ -Curve Loop(142) = {278, 323, 130, -322}; -//+ -Plane Surface(141) = {142}; -//+ -Curve Loop(143) = {279, 324, -129, -323}; -//+ -Plane Surface(142) = {143}; -//+ -Curve Loop(144) = {280, 325, 128, -324}; -//+ -Plane Surface(143) = {144}; -//+ -Curve Loop(145) = {281, 326, -127, -325}; -//+ -Plane Surface(144) = {145}; -//+ -Curve Loop(146) = {282, 327, 126, -326}; -//+ -Plane Surface(145) = {146}; -//+ -//+ -Transfinite Surface {102}; -//+ -Transfinite Surface {103}; -//+ -Transfinite Surface {104}; -//+ -Transfinite Surface {105}; -//+ -Transfinite Surface {106}; -//+ -Transfinite Surface {107}; -//+ -Transfinite Surface {108}; -//+ -Transfinite Surface {109}; -//+ -Transfinite Surface {110}; -//+ -Transfinite Surface {111}; -//+ -Transfinite Surface {112}; -//+ -Transfinite Surface {113}; -//+ -Transfinite Surface {114}; -//+ -Transfinite Surface {115}; -//+ -Transfinite Surface {116}; -//+ -Transfinite Surface {117}; -//+ -Transfinite Surface {118}; -//+ -Transfinite Surface {119}; -//+ -Transfinite Surface {120}; -//+ -Transfinite Surface {121}; -//+ -Transfinite Surface {122}; -//+ -Transfinite Surface {123}; -//+ -Transfinite Surface {124}; -//+ -Transfinite Surface {125}; -//+ -Transfinite Surface {126}; -//+ -Transfinite Surface {127}; -//+ -Transfinite Surface {128}; -//+ -Transfinite Surface {129}; -//+ -Transfinite Surface {130}; -//+ -Transfinite Surface {131}; -//+ -Transfinite Surface {132}; -//+ -Transfinite Surface {133}; -//+ -Transfinite Surface {134}; -//+ -Transfinite Surface {135}; -//+ -Transfinite Surface {136}; -//+ -Transfinite Surface {137}; -//+ -Transfinite Surface {138}; -//+ -Transfinite Surface {139}; -//+ -Transfinite Surface {140}; -//+ -Transfinite Surface {141}; -//+ -Transfinite Surface {142}; -//+ -Transfinite Surface {143}; -//+ -Transfinite Surface {144}; -//+ -Transfinite Surface {145}; -//+ -Transfinite Curve { 238, 239, 240, 241, 242, 243, 244, 245, 250, 251, 252, 253, 254, 255, 256, 257, 262, 263, 264, 265, 266, 267, 268, 269, 274, 275, 276, 277, 278, 279, 280, 281} = n_recombine Using Progression 1; -Transfinite Curve {237,246,249, 258,261,270,273,282} = n_recombine / 4 Using Progression 1; - - -Physical Curve("void",328) = {236, 237, 238, 239, 240, 241, 242, 243, 244, 245, 246, 247, 248, 249, 250, 251, 252, 253, 254, 255, 256, 257, 258, 259, 260, 261, 262, 263, 264, 265, 266, 267, 268, 269, 270, 271, 272, 273, 274, 275, 276, 277, 278, 279, 280, 281, 282, 283}; - -Recombine Surface "*"; -// Define meshing options -//Mesh.Algorithm = 6; // Specify meshing algorithm (e.g., Delaunay) -//Mesh.ElementOrder = 1; // Specify element order -// Generate the mesh -//Mesh 2; \ No newline at end of file diff --git a/examples_hpc/lattice/mesh/make_lattice_mesh.sh b/examples_hpc/lattice/mesh/make_lattice_mesh.sh deleted file mode 100644 index e9bb36e4..00000000 --- a/examples_hpc/lattice/mesh/make_lattice_mesh.sh +++ /dev/null @@ -1,9 +0,0 @@ -#rm lattice_n20.su2 lattice_n20.con -#rm lattice_n40.su2 lattice_n40.con -#rm lattice_n80.su2 lattice_n80.con -rm lattice_n160.su2 lattice_n160.con - -#gmsh lattice_rectangular_n20.geo -2 -format su2 -save_all -o lattice_n20.su2 -#gmsh lattice_rectangular_n40.geo -2 -format su2 -save_all -o lattice_n40.su2 -#gmsh lattice_rectangular_n80.geo -2 -format su2 -save_all -o lattice_n80.su2 -gmsh lattice_rectangular_n160.geo -2 -format su2 -save_all -o lattice_n160.su2 diff --git a/examples_hpc/lattice/run_all.sh b/examples_hpc/lattice/run_all.sh deleted file mode 100644 index b3c665e7..00000000 --- a/examples_hpc/lattice/run_all.sh +++ /dev/null @@ -1,18 +0,0 @@ -#../../build/KiT-RT lattice_S19_n20.cfg -#../../build/KiT-RT lattice_S29_n20.cfg -#../../build/KiT-RT lattice_S41_n20.cfg -#../../build/KiT-RT lattice_S53_n20.cfg -#../../build/KiT-RT lattice_S65_n20.cfg - - -../../build_singularity/KiT-RT lattice_S19_n40.cfg -../../build_singularity/KiT-RT lattice_S29_n40.cfg -../../build_singularity/KiT-RT lattice_S41_n40.cfg -../../build_singularity/KiT-RT lattice_S53_n40.cfg -../../build_singularity/KiT-RT lattice_S65_n40.cfg - -../../build_singularity/KiT-RT lattice_S19_n80.cfg -../../build_singularity/KiT-RT lattice_S29_n80.cfg -../../build_singularity/KiT-RT lattice_S41_n80.cfg -../../build_singularity/KiT-RT lattice_S53_n80.cfg -../../build_singularity/KiT-RT lattice_S65_n80.cfg \ No newline at end of file diff --git a/examples_hpc/sym_hohlraum/sym_hohlraum.cfg b/examples_hpc/sym_hohlraum/sym_hohlraum.cfg index a363f9f1..18092af8 100644 --- a/examples_hpc/sym_hohlraum/sym_hohlraum.cfg +++ b/examples_hpc/sym_hohlraum/sym_hohlraum.cfg @@ -19,18 +19,19 @@ MESH_FILE = mesh/sym_hohlraum_v2_point01.su2 % ---- Problem specifications ---- % PROBLEM = SYMMETRIC_HOHLRAUM -SPATIAL_DIM = 3 +SPATIAL_DIM = 2 % % ---- Design Parameters --- N_SAMPLING_PTS_LINE_GREEN = 100 % % ---- Solver specifications ---- % +HPC_SOLVER = NO CFL_NUMBER = 0.7 TIME_FINAL = 2 SOLVER = SN_SOLVER RECONS_ORDER = 2 -TIME_INTEGRATION_ORDER = 2 +TIME_INTEGRATION_ORDER = 1 % % ---- Boundary Conditions ---- % @@ -38,15 +39,15 @@ BC_NEUMANN = ( void, inflow ) % % ---- Quadrature ---- % -QUAD_TYPE = LEBEDEV -QUAD_ORDER = 9 +QUAD_TYPE = GAUSS_LEGENDRE_TENSORIZED_2D +QUAD_ORDER = 10 % % ----- Output ---- % VOLUME_OUTPUT = (MINIMAL) -VOLUME_OUTPUT_FREQUENCY = 20 +VOLUME_OUTPUT_FREQUENCY = 0 SCREEN_OUTPUT = ( ITER, WALL_TIME, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, TOTAL_PARTICLE_ABSORPTION_CENTER, TOTAL_PARTICLE_ABSORPTION_VERTICAL, TOTAL_PARTICLE_ABSORPTION_HORIZONTAL, PROBE_MOMENT_TIME_TRACE, VAR_ABSORPTION_GREEN) -SCREEN_OUTPUT_FREQUENCY = 20 -HISTORY_OUTPUT = ( ITER, WALL_TIME, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, TOTAL_PARTICLE_ABSORPTION_CENTER, TOTAL_PARTICLE_ABSORPTION_VERTICAL, TOTAL_PARTICLE_ABSORPTION_HORIZONTAL, PROBE_MOMENT_TIME_TRACE, VAR_ABSORPTION_GREEN) +SCREEN_OUTPUT_FREQUENCY = 1 +HISTORY_OUTPUT = ( ITER, WALL_TIME, MASS, RMS_FLUX, VTK_OUTPUT, CSV_OUTPUT, CUR_OUTFLOW, TOTAL_OUTFLOW, MAX_OUTFLOW, TOTAL_PARTICLE_ABSORPTION_CENTER, TOTAL_PARTICLE_ABSORPTION_VERTICAL, TOTAL_PARTICLE_ABSORPTION_HORIZONTAL, PROBE_MOMENT_TIME_TRACE, VAR_ABSORPTION_GREEN, VAR_ABSORPTION_GREEN_LINE) HISTORY_OUTPUT_FREQUENCY = 1 diff --git a/include/common/config.hpp b/include/common/config.hpp index 9c9ca887..52d8c8cb 100644 --- a/include/common/config.hpp +++ b/include/common/config.hpp @@ -48,9 +48,9 @@ class Config // std::vector _1dIntegrationBounds; /*!< @brief Quadrature Order*/ // Mesh - unsigned _nCells; /*!< @brief Number of cells in the mesh */ - unsigned short _dim; /*!< @brief spatial dimensionality of the mesh/test case */ - bool _forcedConnectivityWrite; /*!< @brief If true, the meshconnectivity is always computed and written to .con file */ + unsigned _nCells; /*!< @brief Number of cells in the mesh */ + unsigned short _dim; /*!< @brief spatial dimensionality of the mesh/test case */ + bool _forcedConnectivityWrite = false; /*!< @brief If true, the meshconnectivity is always computed and written to .con file */ // Boundary Conditions /*!< @brief List of all Pairs (marker, BOUNDARY_TYPE), e.g. (farfield,DIRICHLET). @@ -62,6 +62,7 @@ class Config std::vector _MarkerNeumann; /*!< @brief Neumann BC markers. */ // Solver + bool _HPC; /* Triggers better performinc solvers */ double _CFL; /*!< @brief CFL Number for Solver*/ double _tEnd; /*!< @brief Final Time for Simulation */ PROBLEM_NAME _problemName; /*!< @brief Name of predefined Problem */ @@ -304,12 +305,14 @@ class Config bool inline GetForcedConnectivity() const { return _forcedConnectivityWrite; } // Solver Structure + bool inline GetHPC() const { return _HPC; } + double inline GetCFL() const { return _CFL; } bool inline GetCleanFluxMat() const { return _cleanFluxMat; } ENTROPY_NAME inline GetEntropyName() const { return _entropyName; } unsigned short inline GetMaxMomentDegree() const { return _maxMomentDegree; } PROBLEM_NAME inline GetProblemName() const { return _problemName; } - unsigned inline GetReconsOrder() { return _reconsOrder; } + unsigned inline GetSpatialOrder() { return _reconsOrder; } SOLVER_NAME inline GetSolverName() const { return _solverName; } double inline GetTEnd() const { return _tEnd; } bool inline GetSNAllGaussPts() const { return _allGaussPts; } @@ -330,7 +333,7 @@ class Config std::vector inline GetLatticeScatterIndividual() const { return _dsgnScatterIndividual; } unsigned short inline GetNLatticeScatterIndividual() const { return _nDsgnScatterIndividual; } // Hohlraum - unsigned short inline GetNumProbingCellsLineHohlraum() const{ return _nProbingCellsLineGreenHohlraum; } + unsigned short inline GetNumProbingCellsLineHohlraum() const { return _nProbingCellsLineGreenHohlraum; } // Optimizer double inline GetNewtonOptimizerEpsilon() const { return _optimizerEpsilon; } @@ -343,9 +346,9 @@ class Config double inline GetEntropyDynamicAnsatz() const { return _entropyDynamicClosure; } // Neural Closure - unsigned short inline GetModelMK() const{ return _neuralModel; } - unsigned short inline GetNeuralModelGamma() const{ return _neuralGamma; } - bool inline GetEnforceNeuralRotationalSymmetry()const { return _enforceNeuralRotationalSymmetry; } + unsigned short inline GetModelMK() const { return _neuralModel; } + unsigned short inline GetNeuralModelGamma() const { return _neuralGamma; } + bool inline GetEnforceNeuralRotationalSymmetry() const { return _enforceNeuralRotationalSymmetry; } // Boundary Conditions BOUNDARY_TYPE GetBoundaryType( std::string nameMarker ) const; /*!< @brief Get Boundary Type of given marker */ @@ -356,30 +359,30 @@ class Config // Basis name SPHERICAL_BASIS_NAME inline GetSphericalBasisName() const { return _sphericalBasisName; } // Output Structure - std::vector inline GetVolumeOutput() const{ return _volumeOutput; } - unsigned short inline GetNVolumeOutput() const{ return _nVolumeOutput; } - unsigned short inline GetVolumeOutputFrequency() const{ return _volumeOutputFrequency; } + std::vector inline GetVolumeOutput() const { return _volumeOutput; } + unsigned short inline GetNVolumeOutput() const { return _nVolumeOutput; } + unsigned short inline GetVolumeOutputFrequency() const { return _volumeOutputFrequency; } - std::vector inline GetScreenOutput()const { return _screenOutput; } - unsigned short inline GetNScreenOutput() const{ return _nScreenOutput; } - unsigned short inline GetScreenOutputFrequency() const{ return _screenOutputFrequency; } + std::vector inline GetScreenOutput() const { return _screenOutput; } + unsigned short inline GetNScreenOutput() const { return _nScreenOutput; } + unsigned short inline GetScreenOutputFrequency() const { return _screenOutputFrequency; } - std::vector inline GetHistoryOutput() const{ return _historyOutput; } - unsigned short inline GetNHistoryOutput() const{ return _nHistoryOutput; } - unsigned short inline GetHistoryOutputFrequency() const{ return _historyOutputFrequency; } + std::vector inline GetHistoryOutput() const { return _historyOutput; } + unsigned short inline GetNHistoryOutput() const { return _nHistoryOutput; } + unsigned short inline GetHistoryOutputFrequency() const { return _historyOutputFrequency; } // Data generator - bool inline GetDataGeneratorMode()const { return _dataGeneratorMode; } - SAMPLER_NAME inline GetSamplerName() const{ return _sampler; } - unsigned long inline GetTrainingDataSetSize() const{ return _tainingSetSize; } - bool inline GetSizeByDimension() const{ return _sizeByDimension; } - unsigned long inline GetMaxValFirstMoment()const { return _maxValFirstMoment; } // Deprecated - double GetRealizableSetEpsilonU0() const{ return _RealizableSetEpsilonU0; } - double GetRealizableSetEpsilonU1() const{ return _RealizableSetEpsilonU1; } - bool inline GetNormalizedSampling() const{ return _normalizedSampling; } - bool inline GetAlphaSampling() const{ return _alphaSampling; } - bool inline GetUniformSamlping() const{ return _sampleUniform; } - double inline GetAlphaSamplingBound()const { return _alphaBound; } + bool inline GetDataGeneratorMode() const { return _dataGeneratorMode; } + SAMPLER_NAME inline GetSamplerName() const { return _sampler; } + unsigned long inline GetTrainingDataSetSize() const { return _tainingSetSize; } + bool inline GetSizeByDimension() const { return _sizeByDimension; } + unsigned long inline GetMaxValFirstMoment() const { return _maxValFirstMoment; } // Deprecated + double GetRealizableSetEpsilonU0() const { return _RealizableSetEpsilonU0; } + double GetRealizableSetEpsilonU1() const { return _RealizableSetEpsilonU1; } + bool inline GetNormalizedSampling() const { return _normalizedSampling; } + bool inline GetAlphaSampling() const { return _alphaSampling; } + bool inline GetUniformSamlping() const { return _sampleUniform; } + double inline GetAlphaSamplingBound() const { return _alphaBound; } double inline GetMinimalEVBound() const { return _minEVAlphaSampling; } // double inline GetMinimalSamplingVelocity() { return _minSamplingVelocity; } double inline GetMaximalSamplingVelocity() const { return _maxSamplingVelocity; } @@ -387,7 +390,7 @@ class Config double inline GetMaximalSamplingTemperature() const { return _maxSamplingTemperature; } unsigned short inline GetNSamplingTemperatures() const { return _nTemperatures; } bool inline GetIsMomentSolver() const { return _isMomentSolver; } - unsigned short inline GetRKStages() const { return _rungeKuttaStages; } + unsigned short inline GetTemporalOrder() const { return _rungeKuttaStages; } // ---- Setters for option structure // This section is dangerous diff --git a/include/common/linalg.hpp b/include/common/linalg.hpp new file mode 100644 index 00000000..6cf1e877 --- /dev/null +++ b/include/common/linalg.hpp @@ -0,0 +1,35 @@ +/* set of linear algebra functions*/ +#include +#include + +namespace LinAlg { + +inline double dot( std::vector& vec1, std::vector& vec2, unsigned len ) { + double result = 0.0; + for( unsigned i = 0; i < len; i++ ) { + result += vec1[i] * vec2[i]; + } + return result; +} + +inline double l2_norm( std::vector& vec1, unsigned len ) { + double result = 0.0; + for( unsigned i = 0; i < len; i++ ) { + result += vec1[i] * vec1[i]; + } + return sqrt( result ); +} + +inline double mat_vec( std::vector& mat, std::vector& vec1, std::vector& vec2, unsigned n_row, unsigned n_col ) { + /* dim(mat)=n_row x n_col + dim(vec1)=n_row + dim(vec2)=n_col + */ + for( unsigned i = 0; i < n_row; i++ ) { + for( unsigned j = 0; j < n_col; j++ ) { + vec2[j] += mat[i * n_col + j] * vec1[i]; + } + } +} + +} // namespace LinAlg \ No newline at end of file diff --git a/include/problems/lattice.hpp b/include/problems/lattice.hpp index 91554038..e2cd95b6 100644 --- a/include/problems/lattice.hpp +++ b/include/problems/lattice.hpp @@ -21,6 +21,8 @@ class Lattice_SN : public ProblemBase bool IsAbsorption( const Vector& pos ) const; /*!< @return True if pos is in absorption region, False otherwise */ bool IsSource( const Vector& pos ) const; /*!< @return True if pos is in source region, False otherwise */ + bool IsAbsorption( double x, double y ) const; /*!< @return True if pos is in absorption region, False otherwise */ + protected: void SetGhostCells() override; /*!< @brief Sets vector of ghost cells for boundary conditions */ unsigned GetBlockID( const Vector& pos ) const; /*!< @brief Returns checkerboard field id (0-48, row major) of the Lattice test case */ @@ -48,7 +50,7 @@ class Lattice_Moment : public ProblemBase { private: Vector _sigmaS; /*!< @brief Vector of scattering crosssections len: numCells */ - Vector _sigmaT; /*!< @brief Vector of total crosssections. len: numCells*/ + Vector _sigmaT; /*!< @brief Vector of total crosssections. len: numCells*/ Lattice_Moment() = delete; diff --git a/include/solvers/snsolver_hpc.hpp b/include/solvers/snsolver_hpc.hpp new file mode 100644 index 00000000..6521265e --- /dev/null +++ b/include/solvers/snsolver_hpc.hpp @@ -0,0 +1,214 @@ +#ifndef SNSOLVERHPC_H +#define SNSOLVERHPC_H + +// include Matrix, Vector definitions +#include "common/globalconstants.hpp" +#include "common/typedef.hpp" + +// externals +#include "spdlog/spdlog.h" + +// Forward Declarations +class Config; +class Mesh; +class ProblemBase; + +/*! @brief Base class for all solvers. */ +class SNSolverHPC +{ + + private: + double _currTime; /*!< @brief wall-time after current iteration */ + Config* _settings; /*!< @brief config class for global information */ + Mesh* _mesh; + ProblemBase* _problem; + + // Time + unsigned _nIter; /*!< @brief number of time steps, for non CSD, this equals _nEnergies, for _csd, _maxIter =_nEnergies-1*/ + double _dT; /*!< @brief energy/time step size */ + + // Mesh related members, memory optimized + unsigned _nCells; /*!< @brief number of spatial cells */ + unsigned _nSys; /*!< @brief number of equations in the transport system, i.e. num quad pts */ + unsigned _nq; /*!< @brief number of quadrature points */ + unsigned _nDim; + unsigned _nNbr; + unsigned _nNodes; + + std::vector _areas; /*!< @brief surface area of all spatial cells, + dim(_areas) = _NCells */ + std::vector _normals; /*!< @brief edge normals multiplied by edge length, + dim(_normals) = (_NCells,nEdgesPerCell,spatialDim) */ + std::vector _neighbors; /*!< @brief edge neighbor cell ids, dim(_neighbors) = (_NCells,nEdgesPerCell) */ + std::vector _cellMidPoints; /*!< @brief dim _nCells x _dim */ + + std::vector _interfaceMidPoints; /*!< @brief dim: _nCells x _nEdgesPerCell x _dim */ + std::vector _cellBoundaryTypes; /*!< @brief dim: _nCells x _nEdgesPerCell x _dim */ + std::map> _ghostCells; /*!< @brief Vector of ghost cells for boundary conditions. CAN BE MORE EFFICIENT */ + std::vector _relativeInterfaceMidPt; /*!< @brief dim _nCells * _nNbr * _nDim */ + std::vector _relativeCellVertices; /*!< @brief dim _nCells * _nNbr * _nDim */ + + std::map _ghostCellsReflectingY; /*!< map that indicates if a ghostcell has a fixed value or is a mirroring boundary */ + std::vector _quadratureYReflection; /*!< map that gives a Reflection against the y axis for the velocity ordinates */ + + unsigned _temporalOrder; /*!< @brief temporal order (current: 1 & 2) */ + unsigned _spatialOrder; /*!< @brief spatial order (current: 1 & 2) */ + std::vector _solDx; /*!< @brief dim = _nCells x _nSys x _dim*/ + std::vector _limiter; /*!< @brief dim = _nCells x _nSys */ + + // Scattering, absorption and source + std::vector _sigmaS; /*!< @brief dim: _nCells */ + std::vector _sigmaT; /*!< @brief dim: _nCells */ + std::vector _source; /*!< @brief dim: _nCells x _nSys */ + std::vector _scatteringKernel; /*!< @brief dim: _nSys x _nSys */ + + // quadrature related numbers + std::vector _quadPts; /*!< @brief dim: _nSys x _dim*/ + std::vector _quadWeights; /*!< @brief dim: _nSys*/ + + // Solution related members + std::vector _sol; /*!< @brief dim = _nCells x _nSys */ + std::vector _flux; /*!< @brief dim = _nCells x _nSys */ + + // Output related members + std::vector _scalarFlux; /*!< @brief dim = _nCells */ + + // QOIS + double _mass; + double _rmsFlux; + double _curAbsorptionLattice; /*!< @brief Absorption of particles at Lattice checkerboard regions at current time step */ + double _totalAbsorptionLattice; /*!< @brief Absorption of particles at Lattice checkerboard regions integrated until current time step */ + double _curMaxAbsorptionLattice; /*!< @brief Maximum pointwise absorption of particles at Lattice checkerboard regions until current time step */ + double _curScalarOutflow; /*!< @brief Outflow over whole boundary at current time step */ + double _totalScalarOutflow; /*!< @brief Outflow over whole boundary integrated until current time step */ + double _curMaxOrdinateOutflow; /*!< @brief Maximum ordinate-wise ouftlow over boundary over all time steps */ + std::vector _localMaxOrdinateOutflow; /*!< @brief Maximum ordinate-wise ouftlow over boundary over all time steps */ + + // Hohlraum QOIS + double _totalAbsorptionHohlraumCenter; + double _totalAbsorptionHohlraumVertical; + double _totalAbsorptionHohlraumHorizontal; + double _curAbsorptionHohlraumCenter; + double _curAbsorptionHohlraumVertical; + double _curAbsorptionHohlraumHorizontal; + double _varAbsorptionHohlraumGreen; + std::vector _probingCellsHohlraum; /*!< @brief Indices of cells that contain a probing sensor */ + std::vector _probingMoments; /*!< @brief Solution Momnets at the probing cells that contain a probing sensor */ + + unsigned _nProbingCellsLineGreen; /*!< @brief Number of sampling cells that contain a probing sensor for the sliding window */ + std::vector _probingCellsLineGreen; /*!< @brief Indices of cells that contain a probing sensor for the sliding window */ + std::vector _absorptionValsIntegrated; /*!< @brief Avg Absorption value at the sampleing points of lineGreen */ + std::vector _varAbsorptionValsIntegrated; /*!< @brief Var in Avg Absorption value at the sampleing points of lineGreen */ + + // Design parameters + std::vector _cornerUpperLeftGreen; /*!< @brief Coord of corner of the green area (minus thickness/2 of it) relative to the green center */ + std::vector _cornerLowerLeftGreen; /*!< @brief Coord of corner of the green area (minus thickness/2 of it) relative to the green center */ + std::vector + _cornerUpperRightGreen; /*!< @brief Coord of corner of the green area (minus thickness/2 of it) relative to the green center */ + std::vector + _cornerLowerRightGreen; /*!< @brief Coord of corner of the green area (minus thickness/2 of it) relative to the green center */ + + double _widthGreen; /*!< @brief width of the green area */ + double _heightGreen; /*!< @brief height of the green area */ + double _thicknessGreen; /*!< @brief thickness of the green area */ + std::vector _centerGreen; /*!< @brief Center of the Hohlraum */ + + double _redLeftTop; /*!< @brief y coord of the top of the left red area */ + double _redLeftBottom; /*!< @brief y coord of the bottom of the left red area */ + double _redRightTop; /*!< @brief y coord of the top of the right red area */ + double _redRightBottom; /*!< @brief y coord of the bottom of the right red area */ + double _thicknessRedLeft; /*!< @brief thickness of the left red area */ + double _thicknessRedRight; /*!< @brief thickness of the right red area */ + + // Output + std::vector>> _outputFields; /*!< @brief Solver Output: dimensions + (GroupID,FieldID,CellID).*/ + std::vector> _outputFieldNames; /*!< @brief Names of the outputFields: dimensions + (GroupID,FieldID) */ + + std::vector _screenOutputFields; /*!< @brief Solver Output: dimensions (FieldID). */ + std::vector _screenOutputFieldNames; /*!< @brief Names of the outputFields: dimensions (FieldID) */ + + std::vector _historyOutputFields; /*!< @brief Solver Output: dimensions (FieldID). */ + std::vector _historyOutputFieldNames; /*!< @brief Names of the outputFields: dimensions (FieldID) */ + + // ---- Member functions ---- + + // Solver + void FluxOrder1(); + void FluxOrder2(); + + void FVMUpdate(); + + /*! @brief Computes the finite Volume update step for the current iteration + @param idx_iter current (peudo) time iteration */ + void RKUpdate( std::vector& sol0, std::vector& sol_rk ); + + void SetGhostCells(); /*!< @brief Sets vector of ghost cells for + + // Helper + /*! @brief ComputeTimeStep calculates the maximal stable time step using the + cfl number + @param cfl Courant-Friedrichs-Levy condition number */ + double ComputeTimeStep( double cfl ) const; + + // IO + /*! @brief Initializes the output groups and fields of this solver and names + * the fields */ + void PrepareVolumeOutput(); + /*! @brief Function that prepares VTK export and csv export of the current + solver iteration + @param idx_iter current (pseudo) time iteration */ + void WriteVolumeOutput( unsigned idx_iter ); + /*! @brief Save Output solution at given energy (pseudo time) to VTK file. + Write frequency is given by option VOLUME_OUTPUT_FREQUENCY. Always prints + last iteration without iteration affix. + @param idx_iter current (pseudo) time iteration */ + void PrintVolumeOutput( int idx_iter ) const; + /*! @brief Initialized the output fields and their Names for the screenoutput + */ + void PrepareScreenOutput(); + + /*! @brief Function that writes screen and history output fields + @param idx_iter current (pseudo) time iteration */ + void WriteScalarOutput( unsigned idx_iter ); + /*! @brief Prints ScreenOutputFields to Screen and to logger. Write frequency + is given by option SCREEN_OUTPUT_FREQUENCY. Always prints last iteration. + @param idx_iter current (pseudo) time iteration */ + void PrintScreenOutput( unsigned idx_iter ); + /*! @brief Initialized the historyOutputFields and their Names for history + output. Write frequency is given by option HISTORY_OUTPUT_FREQUENCY. Always + prints last iteration. */ + void PrepareHistoryOutput(); + /*! @brief Prints HistoryOutputFields to logger + @param idx_iter current (pseudo) time iteration */ + void PrintHistoryOutput( unsigned idx_iter ); + /*! @brief Pre Solver Screen and Logger Output */ + void DrawPreSolverOutput(); + /*! @brief Post Solver Screen and Logger Output */ + void DrawPostSolverOutput(); + + // Helper + unsigned Idx2D( unsigned idx1, unsigned idx2, unsigned len2 ); + unsigned Idx3D( unsigned idx1, unsigned idx2, unsigned idx3, unsigned len2, unsigned len3 ); + bool IsAbsorptionLattice( double x, double y ) const; + + void SetProbingCellsLineGreen(); + void ComputeQOIsGreenProbingLine(); + std::vector linspace2D( const std::vector& start, const std::vector& end, unsigned num_points ); + + public: + /*! @brief Solver constructor + * @param settings config class that stores all needed config information */ + SNSolverHPC( Config* settings ); + + ~SNSolverHPC(); + + /*! @brief Solve functions runs main iteration loop. Components of the solve + * loop are pure and subclassed by the child solvers. */ + void Solve(); + + /*! @brief Save Output solution to VTK file */ + void PrintVolumeOutput() const; // Only for debugging purposes. +}; +#endif // SNSOLVERHPC_H diff --git a/include/solvers/solverbase.hpp b/include/solvers/solverbase.hpp index 47d53861..8b6e7551 100644 --- a/include/solvers/solverbase.hpp +++ b/include/solvers/solverbase.hpp @@ -75,7 +75,7 @@ class SolverBase // _nCells x nSys. nSys depents on the solver (either moments or quadpts) VectorVector _solDx; /*!< @brief solution gradient in x direction*/ VectorVector _solDy; /*!< @brief solution gradient in y direction*/ - VectorVector _limiter; /*! < @brief slope limiter at cell and system pt */ + VectorVector _limiter; /*!< @brief slope limiter at cell and system pt */ // Output related members std::vector>> _outputFields; /*!< @brief Solver Output: dimensions @@ -108,7 +108,7 @@ class SolverBase virtual void FVMUpdate( unsigned idx_iter ) = 0; /*! @brief Computes the finite Volume update step for the current iteration @param idx_iter current (peudo) time iteration */ - virtual void RKUpdate( VectorVector sol0, VectorVector sol_rk ); + virtual void RKUpdate( VectorVector& sol0, VectorVector& sol_rk ); // Helper /*! @brief ComputeTimeStep calculates the maximal stable time step using the diff --git a/include/toolboxes/reconstructor.hpp b/include/toolboxes/reconstructor.hpp index 4f88d629..ee126369 100644 --- a/include/toolboxes/reconstructor.hpp +++ b/include/toolboxes/reconstructor.hpp @@ -27,7 +27,7 @@ class Reconstructor virtual ~Reconstructor() {} static Reconstructor* Create( Config* settings ); - unsigned inline GetReconsOrder() { return _reconsOrder; } + unsigned inline GetSpatialOrder() { return _reconsOrder; } /*! Method 1: structured developing * @brief Slope of angular flux psi inside a given cell diff --git a/src/common/config.cpp b/src/common/config.cpp index 73c04200..50124351 100644 --- a/src/common/config.cpp +++ b/src/common/config.cpp @@ -231,6 +231,9 @@ void Config::SetConfigOptions() { AddUnsignedShortOption( "QUAD_ORDER", _quadOrder, 1 ); // Solver related options + /*! @brief HPC_SOLVER \n DESCRIPTION: If true, the SN Solver uses hpc implementation. \n DEFAULT false \ingroup Config */ + AddBoolOption( "HPC_SOLVER", _HPC, false ); + /*! @brief MAX_MOMENT_ORDER \n: DESCRIPTON: Specifies the maximal order of Moments for PN and SN Solver */ AddUnsignedShortOption( "MAX_MOMENT_SOLVER", _maxMomentDegree, 1 ); /*! @brief CFL \n DESCRIPTION: CFL number \n DEFAULT 1.0 @ingroup Config.*/ @@ -252,7 +255,7 @@ void Config::SetConfigOptions() { /*! @brief Runge Kutta Staes \n DESCRIPTION: Sets number of Runge Kutta Stages for time integration \n DEFAULT 1 \ingroup Config */ AddUnsignedShortOption( "TIME_INTEGRATION_ORDER", _rungeKuttaStages, 1 ); - // Problem Relateed Options + // Problem Related Options /*! @brief MaterialDir \n DESCRIPTION: Relative Path to the data directory (used in the ICRU database class), starting from the directory of the * cfg file . \n DEFAULT "../data/material/" \ingroup Config */ AddStringOption( "DATA_DIR", _dataDir, string( "../data/" ) ); @@ -578,7 +581,7 @@ void Config::SetPostprocessing() { ErrorMessages::Error( "CSD_MN_SOLVER only works with Spherical Harmonics currently.", CURRENT_FUNCTION ); } - if( GetReconsOrder() > 2 ) { + if( GetSpatialOrder() > 2 ) { ErrorMessages::Error( "Solvers only work with 1st and 2nd order spatial fluxes.", CURRENT_FUNCTION ); } @@ -626,10 +629,10 @@ void Config::SetPostprocessing() { for( unsigned short idx_volOutput = 0; idx_volOutput < _nVolumeOutput; idx_volOutput++ ) { switch( _solverName ) { case SN_SOLVER: - if( _problemName == PROBLEM_SymmetricHohlraum ) + if( _problemName == PROBLEM_Linesource ) supportedGroups = { MINIMAL, ANALYTIC }; else - supportedGroups = { MINIMAL, ANALYTIC }; + supportedGroups = { MINIMAL }; if( supportedGroups.end() == std::find( supportedGroups.begin(), supportedGroups.end(), _volumeOutput[idx_volOutput] ) ) { ErrorMessages::Error( "SN_SOLVER only supports volume output MINIMAL and ANALYTIC.\nPlease check your .cfg file.", @@ -988,6 +991,16 @@ void Config::SetPostprocessing() { for( unsigned i = 0; i < 6; i++ ) _historyOutput.push_back( PROBE_MOMENT_TIME_TRACE ); } } + + if( _problemName == PROBLEM_SymmetricHohlraum ) { + std::vector::iterator it; + it = find( _historyOutput.begin(), _historyOutput.end(), VAR_ABSORPTION_GREEN_LINE ); + if( it != _historyOutput.end() ) { + _historyOutput.erase( it ); + _nHistoryOutput += _nProbingCellsLineGreenHohlraum - 1; // extend the screen output by the number of probing points + for( unsigned i = 0; i < _nProbingCellsLineGreenHohlraum; i++ ) _historyOutput.push_back( VAR_ABSORPTION_GREEN_LINE ); + } + } } // Mesh postprocessing diff --git a/src/common/io.cpp b/src/common/io.cpp index 61224c62..1debf80d 100644 --- a/src/common/io.cpp +++ b/src/common/io.cpp @@ -336,6 +336,7 @@ void LoadConnectivityFromFile( const std::string inputFile, unsigned correctedNodesPerCell = nNodesPerCell; if( count < nNodesPerCell + nNodesPerCell * nDim * 2 + 1 ) { correctedNodesPerCell = nNodesPerCell - 1; + ErrorMessages::Error( "Error opening connectivity file.", CURRENT_FUNCTION ); } std::istringstream iss( line ); diff --git a/src/common/mesh.cpp b/src/common/mesh.cpp index 6f82b537..7cf3f8e8 100644 --- a/src/common/mesh.cpp +++ b/src/common/mesh.cpp @@ -191,18 +191,20 @@ void Mesh::ComputeConnectivity() { unsigned IDj = neighborsFlat[i]; if( IDi == IDj ) continue; // avoid self assignment if( IDj == _ghostCellID ) { // cell is boundary cell - if( std::find( _cellNeighbors[IDi].begin(), _cellNeighbors[IDi].end(), _ghostCellID ) == _cellNeighbors[IDi].end() ) { - _cellNeighbors[IDi].push_back( _ghostCellID ); - _cellNormals[IDi].push_back( normalsFlat[i] ); - _cellInterfaceMidPoints[IDi].push_back( interfaceMidFlat[i] ); - } + // if( std::find( _cellNeighbors[IDi].begin(), _cellNeighbors[IDi].end(), _ghostCellID ) == _cellNeighbors[IDi].end() ) { + _cellNeighbors[IDi].push_back( _ghostCellID ); + _cellNormals[IDi].push_back( normalsFlat[i] ); + _cellInterfaceMidPoints[IDi].push_back( interfaceMidFlat[i] ); + // temp++; + // } + // std::cout << temp << "\n"; } else { // normal cell neighbor - if( std::find( _cellNeighbors[IDi].begin(), _cellNeighbors[IDi].end(), IDj ) == _cellNeighbors[IDi].end() ) { - _cellNeighbors[IDi].push_back( neighborsFlat[i] ); - _cellNormals[IDi].push_back( normalsFlat[i] ); - _cellInterfaceMidPoints[IDi].push_back( interfaceMidFlat[i] ); - } + // if( std::find( _cellNeighbors[IDi].begin(), _cellNeighbors[IDi].end(), IDj ) == _cellNeighbors[IDi].end() ) { + _cellNeighbors[IDi].push_back( neighborsFlat[i] ); + _cellNormals[IDi].push_back( normalsFlat[i] ); + _cellInterfaceMidPoints[IDi].push_back( interfaceMidFlat[i] ); + //} } } @@ -361,15 +363,15 @@ void Mesh::ComputeLimiter( // Compute value at interface midpoint, called gaussPt double gaussPt = 0.0; // gauss point is at cell vertex - gaussPt = 0.5 * ( solDx[idx_cell][idx_sys] * ( _nodes[_cells[idx_cell][idx_nbr]][0] - _cellMidPoints[idx_cell][0] ) + - solDy[idx_cell][idx_sys] * ( _nodes[_cells[idx_cell][idx_nbr]][1] - _cellMidPoints[idx_cell][1] ) ); + gaussPt = ( solDx[idx_cell][idx_sys] * ( _nodes[_cells[idx_cell][idx_nbr]][0] - _cellMidPoints[idx_cell][0] ) + + solDy[idx_cell][idx_sys] * ( _nodes[_cells[idx_cell][idx_nbr]][1] - _cellMidPoints[idx_cell][1] ) ); // Compute limiter input if( std::abs( gaussPt ) > eps ) { if( gaussPt > 0.0 ) { r = ( maxSol - sol[idx_cell][idx_sys] ) / gaussPt; } - else if( gaussPt < 0.0 ) { + else { r = ( minSol - sol[idx_cell][idx_sys] ) / gaussPt; } } diff --git a/src/fluxes/upwindflux.cpp b/src/fluxes/upwindflux.cpp index 8f3e9041..2caac866 100644 --- a/src/fluxes/upwindflux.cpp +++ b/src/fluxes/upwindflux.cpp @@ -14,14 +14,13 @@ double UpwindFlux::Flux( const Vector& Omega, double psiL, double psiR, const Ve } void UpwindFlux::Flux( const VectorVector& quadPts, const Vector& psiL, const Vector& psiR, Vector& flux, const Vector& n, unsigned n_sys ) { - double inner; // Only use x and y axis in 2d case, minus because y axis is flipped for some reason + double inner; // Only use x and y axis in 2d case, minus because y axis is flipped for some reason - for (unsigned idx_q = 0; idx_q < n_sys; idx_q++){ - inner = quadPts[idx_q][0] * n[0] + quadPts[idx_q][1] * n[1]; + for( unsigned idx_q = 0; idx_q < n_sys; idx_q++ ) { + inner = quadPts[idx_q][0] * n[0] + quadPts[idx_q][1] * n[1]; - flux[idx_q] += (inner > 0) ? inner * psiL[idx_q] : inner * psiR[idx_q]; - - } + flux[idx_q] += ( inner > 0 ) ? inner * psiL[idx_q] : inner * psiR[idx_q]; + } } double UpwindFlux::FluxXZ( const Vector& Omega, double psiL, double psiR, const Vector& n ) const { diff --git a/src/main.cpp b/src/main.cpp index 4390e03a..57e1f389 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -10,7 +10,7 @@ #include "common/config.hpp" #include "common/io.hpp" #include "solvers/solverbase.hpp" - +#include "solvers/snsolver_hpc.hpp" #include "datagenerator/datageneratorbase.hpp" #ifdef BUILD_GUI @@ -45,11 +45,20 @@ int main( int argc, char** argv ) { } else { // Build solver - SolverBase* solver = SolverBase::Create( config ); - // Run solver and export - solver->Solve(); - solver->PrintVolumeOutput(); - delete solver; + if (config->GetHPC()){ + SNSolverHPC* solver = new SNSolverHPC( config ); + // Run solver and export + solver->Solve(); + solver->PrintVolumeOutput(); + delete solver; + } + else{ + SolverBase* solver = SolverBase::Create( config ); + // Run solver and export + solver->Solve(); + solver->PrintVolumeOutput(); + delete solver; + } } delete config; diff --git a/src/problems/halflattice.cpp b/src/problems/halflattice.cpp index c4279d18..ba994093 100644 --- a/src/problems/halflattice.cpp +++ b/src/problems/halflattice.cpp @@ -71,7 +71,7 @@ std::vector HalfLattice_SN::GetExternalSource( const Vector& /*ene } VectorVector HalfLattice_SN::SetupIC() { - VectorVector psi( _mesh->GetNumCells(), Vector( _settings->GetNQuadPoints(), 1e-10 ) ); + VectorVector psi( _mesh->GetNumCells(), Vector( _settings->GetNQuadPoints(), 1e-15 ) ); return psi; } @@ -198,6 +198,7 @@ const Vector& HalfLattice_SN::GetGhostCellValue( int idx_cell, const Vector& cel double HalfLattice_SN::GetCurAbsorptionLattice() { return _curAbsorptionLattice; } double HalfLattice_SN::GetTotalAbsorptionLattice() { return _totalAbsorptionLattice; } double HalfLattice_SN::GetMaxAbsorptionLattice() { return _curMaxAbsorptionLattice; } + // QOI setter void HalfLattice_SN::ComputeTotalAbsorptionLattice( double dT ) { _totalAbsorptionLattice += _curAbsorptionLattice * dT; } diff --git a/src/problems/lattice.cpp b/src/problems/lattice.cpp index 769e1b3a..4c258c5b 100644 --- a/src/problems/lattice.cpp +++ b/src/problems/lattice.cpp @@ -71,7 +71,7 @@ std::vector Lattice_SN::GetExternalSource( const Vector& /*energie } VectorVector Lattice_SN::SetupIC() { - VectorVector psi( _mesh->GetNumCells(), Vector( _settings->GetNQuadPoints(), 1e-10 ) ); + VectorVector psi( _mesh->GetNumCells(), Vector( _settings->GetNQuadPoints(), 1e-15 ) ); return psi; } @@ -91,6 +91,22 @@ bool Lattice_SN::IsAbsorption( const Vector& pos ) const { return false; } +bool Lattice_SN::IsAbsorption( double x, double y ) const { + // Check whether pos is inside absorbing squares + double xy_corrector = -3.5; + std::vector lbounds{ 1 + xy_corrector, 2 + xy_corrector, 3 + xy_corrector, 4 + xy_corrector, 5 + xy_corrector }; + std::vector ubounds{ 2 + xy_corrector, 3 + xy_corrector, 4 + xy_corrector, 5 + xy_corrector, 6 + xy_corrector }; + for( unsigned k = 0; k < lbounds.size(); ++k ) { + for( unsigned l = 0; l < lbounds.size(); ++l ) { + if( ( l + k ) % 2 == 1 || ( k == 2 && l == 2 ) || ( k == 2 && l == 4 ) ) continue; + if( x >= lbounds[k] && x <= ubounds[k] && y >= lbounds[l] && y <= ubounds[l] ) { + return true; + } + } + } + return false; +} + unsigned Lattice_SN::GetBlockID( const Vector& pos ) const { double xy_corrector = 3.5; int block_x = int( pos[0] + xy_corrector ); @@ -119,7 +135,7 @@ void Lattice_SN::SetGhostCells() { for( unsigned idx_cell = 0; idx_cell < _mesh->GetNumCells(); idx_cell++ ) { if( cellBoundaries[idx_cell] == BOUNDARY_TYPE::NEUMANN || cellBoundaries[idx_cell] == BOUNDARY_TYPE::DIRICHLET ) { - ghostCellMap.insert( { idx_cell, void_ghostcell } ); + ghostCellMap[idx_cell] = void_ghostcell; } } _ghostCells = ghostCellMap; diff --git a/src/problems/linesource.cpp b/src/problems/linesource.cpp index 94d70b61..8fc22ac9 100644 --- a/src/problems/linesource.cpp +++ b/src/problems/linesource.cpp @@ -129,7 +129,7 @@ VectorVector LineSource_SN::SetupIC() { for( unsigned j = 0; j < cellMids.size(); ++j ) { double x = cellMids[j][0]; double y = cellMids[j][1]; - psi[j] = Vector( _settings->GetNQuadPoints(), 1.0 / ( 4.0 * M_PI * t ) * std::exp( -( x * x + y * y ) / ( 4 * t ) ) ) / ( 4 * M_PI ); + psi[j] = Vector( _settings->GetNQuadPoints(), 10.0 / ( 4.0 * M_PI * t ) * std::exp( -( x * x + y * y ) / ( 4 * t ) ) ) / ( 4 * M_PI ); } return psi; } diff --git a/src/quadratures/qgausslegendretensorized2D.cpp b/src/quadratures/qgausslegendretensorized2D.cpp index 5f26bfd6..b1fc62cf 100644 --- a/src/quadratures/qgausslegendretensorized2D.cpp +++ b/src/quadratures/qgausslegendretensorized2D.cpp @@ -58,7 +58,7 @@ void QGaussLegendreTensorized2D::SetPointsAndWeights() { } unsigned range = std::floor( _order / 2.0 ); // comment (steffen): Only half of the points, due to projection - double normalizationFactor = .5; + double normalizationFactor = 2; // resize points and weights _pointsKarth.resize( _nq ); diff --git a/src/solvers/snsolver.cpp b/src/solvers/snsolver.cpp index a5c5e801..527e48cb 100644 --- a/src/solvers/snsolver.cpp +++ b/src/solvers/snsolver.cpp @@ -112,28 +112,28 @@ void SNSolver::FluxUpdatePseudo2D() { // Loop over all spatial cells #pragma omp parallel for for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) { - double solL; - double solR; // Dirichlet cells stay at IC, farfield assumption if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::DIRICHLET ) continue; - // Reset temporary variable - _solNew[idx_cell]*= 0.0; //blaze op + // Reset temporary variable + _solNew[idx_cell] *= 0.0; // blaze op - // Loop over all neighbor cells (edges) of cell j and compute numerical - // fluxes - for( unsigned idx_nbr = 0; idx_nbr < _neighbors[idx_cell].size(); ++idx_nbr ) { - // store flux contribution on psiNew_sigmaS to save memory - if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::NEUMANN && _neighbors[idx_cell][idx_nbr] == _nCells ) { - _g->Flux( _quadPoints, _sol[idx_cell], _problem->GetGhostCellValue( idx_cell, _sol[idx_cell] ), - _solNew[idx_cell], - _normals[idx_cell][idx_nbr], _nq ); - } - else { - // first order solver - _g->Flux( _quadPoints, _sol[idx_cell], _sol[ _neighbors[idx_cell][idx_nbr]], _solNew[idx_cell],_normals[idx_cell][idx_nbr], _nq ); - } + // Loop over all neighbor cells (edges) of cell j and compute numerical + // fluxes + for( unsigned idx_nbr = 0; idx_nbr < _neighbors[idx_cell].size(); ++idx_nbr ) { + // store flux contribution on psiNew_sigmaS to save memory + if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::NEUMANN && _neighbors[idx_cell][idx_nbr] == _nCells ) { + _g->Flux( _quadPoints, + _sol[idx_cell], + _problem->GetGhostCellValue( idx_cell, _sol[idx_cell] ), + _solNew[idx_cell], + _normals[idx_cell][idx_nbr], + _nq ); + } + else { + // first order solver + _g->Flux( _quadPoints, _sol[idx_cell], _sol[_neighbors[idx_cell][idx_nbr]], _solNew[idx_cell], _normals[idx_cell][idx_nbr], _nq ); } - + } } } else if( _reconsOrder == 2 ) { @@ -145,30 +145,30 @@ void SNSolver::FluxUpdatePseudo2D() { // Dirichlet cells stay at IC, farfield assumption if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::DIRICHLET ) continue; // Loop over all ordinates - _solNew[idx_cell] *= 0.0; // blaze operation - // Loop over all neighbor cells (edges) of cell j and compute numerical - // fluxes + _solNew[idx_cell] *= 0.0; // blaze operation + // Loop over all neighbor cells (edges) of cell j and compute numerical + // fluxes for( unsigned idx_nbr = 0; idx_nbr < _neighbors[idx_cell].size(); ++idx_nbr ) { // store flux contribution on psiNew_sigmaS to save memory if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::NEUMANN && _neighbors[idx_cell][idx_nbr] == _nCells ) { - _g->Flux( _quadPoints, _sol[idx_cell], _problem->GetGhostCellValue( idx_cell, _sol[idx_cell] ), - _solNew[idx_cell], - _normals[idx_cell][idx_nbr], _nq ); + _g->Flux( _quadPoints, + _sol[idx_cell], + _problem->GetGhostCellValue( idx_cell, _sol[idx_cell] ), + _solNew[idx_cell], + _normals[idx_cell][idx_nbr], + _nq ); } else { unsigned int nbr_glob = _neighbors[idx_cell][idx_nbr]; // global idx of neighbor cell // left status of interface - solL = _sol[idx_cell] + - _limiter[idx_cell] * - (_solDx[idx_cell] * ( _interfaceMidPoints[idx_cell][idx_nbr][0] - _cellMidPoints[idx_cell][0] ) + - _solDy[idx_cell] * ( _interfaceMidPoints[idx_cell][idx_nbr][1] - _cellMidPoints[idx_cell][1] ) ); - + solL = _sol[idx_cell] + + _limiter[idx_cell] * ( _solDx[idx_cell] * ( _interfaceMidPoints[idx_cell][idx_nbr][0] - _cellMidPoints[idx_cell][0] ) + + _solDy[idx_cell] * ( _interfaceMidPoints[idx_cell][idx_nbr][1] - _cellMidPoints[idx_cell][1] ) ); solR = _sol[nbr_glob] + - _limiter[nbr_glob] * - ( _solDx[nbr_glob]* ( _interfaceMidPoints[idx_cell][idx_nbr][0] - _cellMidPoints[nbr_glob][0] ) + - _solDy[nbr_glob] * ( _interfaceMidPoints[idx_cell][idx_nbr][1] - _cellMidPoints[nbr_glob][1] ) ); + _limiter[nbr_glob] * ( _solDx[nbr_glob] * ( _interfaceMidPoints[idx_cell][idx_nbr][0] - _cellMidPoints[nbr_glob][0] ) + + _solDy[nbr_glob] * ( _interfaceMidPoints[idx_cell][idx_nbr][1] - _cellMidPoints[nbr_glob][1] ) ); // flux evaluation _g->Flux( _quadPoints, solL, solR, _solNew[idx_cell], _normals[idx_cell][idx_nbr], _nq ); @@ -185,19 +185,19 @@ void SNSolver::FVMUpdate( unsigned idx_iter ) { // Dirichlet cells stay at IC, farfield assumption if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::DIRICHLET ) continue; // loop over all ordinates - //for( unsigned idx_quad = 0; idx_quad < _nq; ++idx_quad ) { - // time update angular flux with numerical flux and total scattering cross - // section - _solNew[idx_cell] = _sol[idx_cell] - ( _dT / _areas[idx_cell] ) * _solNew[idx_cell] - - _dT * _sigmaT[0][idx_cell] * _sol[idx_cell]; + // for( unsigned idx_quad = 0; idx_quad < _nq; ++idx_quad ) { + // time update angular flux with numerical flux and total scattering cross + // section + + _solNew[idx_cell] = _sol[idx_cell] - ( _dT / _areas[idx_cell] ) * _solNew[idx_cell] - _dT * _sigmaT[0][idx_cell] * _sol[idx_cell]; //} // compute scattering effects _solNew[idx_cell] += _dT * _sigmaS[0][idx_cell] * _scatteringKernel * _sol[idx_cell]; // multiply scattering matrix with psi // Source Term - //if( _Q[0][idx_cell].size() == 1u ) // isotropic source + // if( _Q[0][idx_cell].size() == 1u ) // isotropic source _solNew[idx_cell] += _dT * _Q[0][idx_cell][0]; - //else - // _solNew[idx_cell] += _dT * _Q[0][idx_cell]; + // else + // _solNew[idx_cell] += _dT * _Q[0][idx_cell]; } } else { @@ -206,11 +206,10 @@ void SNSolver::FVMUpdate( unsigned idx_iter ) { // Dirichlet cells stay at IC, farfield assumption if( _boundaryCells[idx_cell] == BOUNDARY_TYPE::DIRICHLET ) continue; // loop over all ordinates - //for( unsigned idx_quad = 0; idx_quad < _nq; ++idx_quad ) { - // time update angular flux with numerical flux and total scattering cross - // section - _solNew[idx_cell] = _sol[idx_cell] - ( _dT / _areas[idx_cell] ) * _solNew[idx_cell] - - _dT * _sigmaT[idx_iter][idx_cell] * _sol[idx_cell]; + // for( unsigned idx_quad = 0; idx_quad < _nq; ++idx_quad ) { + // time update angular flux with numerical flux and total scattering cross + // section + _solNew[idx_cell] = _sol[idx_cell] - ( _dT / _areas[idx_cell] ) * _solNew[idx_cell] - _dT * _sigmaT[idx_iter][idx_cell] * _sol[idx_cell]; // } // compute scattering effects _solNew[idx_cell] += _dT * _sigmaS[idx_iter][idx_cell] * _scatteringKernel * _sol[idx_cell]; // multiply scattering matrix with psi diff --git a/src/solvers/snsolver_hpc.cpp b/src/solvers/snsolver_hpc.cpp new file mode 100644 index 00000000..ef9a7f63 --- /dev/null +++ b/src/solvers/snsolver_hpc.cpp @@ -0,0 +1,1262 @@ +#include "solvers/snsolver_hpc.hpp" +#include "common/config.hpp" +#include "common/io.hpp" +#include "common/mesh.hpp" +#include "kernels/scatteringkernelbase.hpp" +#include "problems/problembase.hpp" +#include "quadratures/quadraturebase.hpp" +#include "toolboxes/textprocessingtoolbox.hpp" + +SNSolverHPC::SNSolverHPC( Config* settings ) { + _settings = settings; + _currTime = 0.0; + + // Create Mesh + _mesh = LoadSU2MeshFromFile( settings ); + _settings->SetNCells( _mesh->GetNumCells() ); + auto quad = QuadratureBase::Create( settings ); + _settings->SetNQuadPoints( quad->GetNq() ); + + _problem = ProblemBase::Create( _settings, _mesh, quad ); + + _nCells = _mesh->GetNumCells(); + _nNbr = _mesh->GetNumNodesPerCell(); + _nDim = _mesh->GetDim(); + _nNodes = _mesh->GetNumNodes(); + + _nq = quad->GetNq(); + _nSys = _nq; + + _spatialOrder = _settings->GetSpatialOrder(); + _temporalOrder = _settings->GetTemporalOrder(); + + _areas = std::vector( _nCells ); + _normals = std::vector( _nCells * _nNbr * _nDim ); + _neighbors = std::vector( _nCells * _nNbr ); + _cellMidPoints = std::vector( _nCells * _nDim ); + _interfaceMidPoints = std::vector( _nCells * _nNbr * _nDim ); + _cellBoundaryTypes = std::vector( _nCells ); + _relativeInterfaceMidPt = std::vector( _nCells * _nNbr * _nDim ); + _relativeCellVertices = std::vector( _nCells * _nNbr * _nDim ); + + // Slope + _solDx = std::vector( _nCells * _nSys * _nDim ); + _limiter = std::vector( _nCells * _nSys ); + + // Physics + _sigmaS = std::vector( _nCells ); + _sigmaT = std::vector( _nCells ); + _source = std::vector( _nCells * _nSys ); + _scatteringKernel = std::vector( _nSys * _nSys ); + + // Quadrature + _quadPts = std::vector( _nSys * _nDim ); + _quadWeights = std::vector( _nSys ); + _scatteringKernel = std::vector( _nSys * _nSys ); + _quadratureYReflection = std::vector( _nSys ); + + // Solution + _sol = std::vector( _nCells * _nSys ); + _flux = std::vector( _nCells * _nSys ); + + _scalarFlux = std::vector( _nCells ); + _localMaxOrdinateOutflow = std::vector( _nCells ); + + auto areas = _mesh->GetCellAreas(); + auto neighbors = _mesh->GetNeighbours(); + auto normals = _mesh->GetNormals(); + auto cellMidPts = _mesh->GetCellMidPoints(); + auto interfaceMidPts = _mesh->GetInterfaceMidPoints(); + auto boundaryTypes = _mesh->GetBoundaryTypes(); + auto nodes = _mesh->GetNodes(); + auto cells = _mesh->GetCells(); + + for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) { + _areas[idx_cell] = areas[idx_cell]; + } + + _dT = ComputeTimeStep( _settings->GetCFL() ); + _nIter = unsigned( _settings->GetTEnd() / _dT ) + 1; + + auto quadPoints = quad->GetPoints(); + auto quadWeights = quad->GetWeights(); + + auto initialCondition = _problem->SetupIC(); + auto sigmaT = _problem->GetTotalXS( Vector( _nIter, 0.0 ) ); + auto sigmaS = _problem->GetScatteringXS( Vector( _nIter, 0.0 ) ); + auto source = _problem->GetExternalSource( Vector( _nIter, 0.0 ) ); + + // Copy to everything to solver + _mass = 0; + + ScatteringKernel* k = ScatteringKernel::CreateScatteringKernel( settings->GetKernelName(), quad ); + auto scatteringKernel = k->GetScatteringKernel(); + for( unsigned idx_sys = 0; idx_sys < _nSys; idx_sys++ ) { + for( unsigned idx_dim = 0; idx_dim < _nDim; idx_dim++ ) { + _quadPts[Idx2D( idx_sys, idx_dim, _nDim )] = quadPoints[idx_sys][idx_dim]; + } + _quadWeights[idx_sys] = quadWeights[idx_sys]; + + for( unsigned idx_sys2 = 0; idx_sys2 < _nSys; idx_sys2++ ) { + _scatteringKernel[Idx2D( idx_sys, idx_sys2, _nSys )] = scatteringKernel( idx_sys, idx_sys2 ); + } + } + + // #pragma omp parallel for + for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) { + _cellBoundaryTypes[idx_cell] = boundaryTypes[idx_cell]; + + for( unsigned idx_dim = 0; idx_dim < _nDim; idx_dim++ ) { + _cellMidPoints[Idx2D( idx_cell, idx_dim, _nDim )] = cellMidPts[idx_cell][idx_dim]; + } + + for( unsigned idx_nbr = 0; idx_nbr < _nNbr; idx_nbr++ ) { + + _neighbors[Idx2D( idx_cell, idx_nbr, _nNbr )] = neighbors[idx_cell][idx_nbr]; + + for( unsigned idx_dim = 0; idx_dim < _nDim; idx_dim++ ) { + _normals[Idx3D( idx_cell, idx_nbr, idx_dim, _nNbr, _nDim )] = normals[idx_cell][idx_nbr][idx_dim]; + _interfaceMidPoints[Idx3D( idx_cell, idx_nbr, idx_dim, _nNbr, _nDim )] = interfaceMidPts[idx_cell][idx_nbr][idx_dim]; + _relativeInterfaceMidPt[Idx3D( idx_cell, idx_nbr, idx_dim, _nNbr, _nDim )] = + _interfaceMidPoints[Idx3D( idx_cell, idx_nbr, idx_dim, _nNbr, _nDim )] - _cellMidPoints[Idx2D( idx_cell, idx_dim, _nDim )]; + _relativeCellVertices[Idx3D( idx_cell, idx_nbr, idx_dim, _nNbr, _nDim )] = + nodes[cells[idx_cell][idx_nbr]][idx_dim] - cellMidPts[idx_cell][idx_dim]; + } + } + + _sigmaS[idx_cell] = sigmaS[0][idx_cell]; + _sigmaT[idx_cell] = sigmaT[0][idx_cell]; + _scalarFlux[idx_cell] = 0; + + for( unsigned idx_sys = 0; idx_sys < _nSys; idx_sys++ ) { + _source[Idx2D( idx_cell, idx_sys, _nSys )] = source[0][idx_cell][0]; // CAREFUL HERE hardcoded to isotropic source + _sol[Idx2D( idx_cell, idx_sys, _nSys )] = initialCondition[idx_cell][idx_sys]; + + _scalarFlux[idx_cell] += _sol[Idx2D( idx_cell, idx_sys, _nSys )] * _quadWeights[idx_sys]; + } + _mass += _scalarFlux[idx_cell] * _areas[idx_cell]; + } + + SetGhostCells(); // ONLY FOR LATTICE AND HALF LATTICE + + PrepareScreenOutput(); // Screen Output + PrepareHistoryOutput(); // History Output + + delete quad; + delete k; + + // Initialiye QOIS + _mass = 0; + _rmsFlux = 0; + + // Lattice + _curAbsorptionLattice = 0; + _totalAbsorptionLattice = 0; + _curMaxAbsorptionLattice = 0; + _curScalarOutflow = 0; + _totalScalarOutflow = 0; + _curMaxOrdinateOutflow = 0; + + // Hohlraum + _totalAbsorptionHohlraumCenter = 0; + _totalAbsorptionHohlraumVertical = 0; + _totalAbsorptionHohlraumHorizontal = 0; + _curAbsorptionHohlraumCenter = 0; + _curAbsorptionHohlraumVertical = 0; + _curAbsorptionHohlraumHorizontal = 0; + _varAbsorptionHohlraumGreen = 0; + + _probingCellsHohlraum = { + _mesh->GetCellOfKoordinate( -0.4, 0. ), + _mesh->GetCellOfKoordinate( 0.4, 0. ), + _mesh->GetCellOfKoordinate( 0., -0.6 ), + _mesh->GetCellOfKoordinate( 0., 0.6 ), + }; + _probingMoments = std::vector( 12, 0. ); + + // Red + _redLeftTop = 0.4; + _redLeftBottom = -0.4; + _redRightTop = 0.4; + _redRightBottom = -0.4; + _thicknessRedLeft = 0.05; + _thicknessRedRight = 0.05; + // Green + _widthGreen = 0.4; + _heightGreen = 0.8; + _thicknessGreen = 0.05; + _centerGreen = { 0.0, 0.0 }; + + _cornerUpperLeftGreen = { -0.2 + _thicknessGreen / 2.0, 0.4 - _thicknessGreen / 2.0 }; + _cornerLowerLeftGreen = { -0.2 + _thicknessGreen / 2.0, -0.4 + _thicknessGreen / 2.0 }; + _cornerUpperRightGreen = { 0.2 - _thicknessGreen / 2.0, 0.4 - _thicknessGreen / 2.0 }; + _cornerLowerRightGreen = { 0.2 - _thicknessGreen / 2.0, -0.4 + _thicknessGreen / 2.0 }; + + std::cout << "Solver initialized!" << std::endl; + + _nProbingCellsLineGreen = _settings->GetNumProbingCellsLineHohlraum(); + SetProbingCellsLineGreen(); + std::cout << "Probing cells initialized!" << std::endl; + + _absorptionValsIntegrated = std::vector( _nProbingCellsLineGreen, 0.0 ); + _varAbsorptionValsIntegrated = std::vector( _nProbingCellsLineGreen, 0.0 ); +} + +SNSolverHPC::~SNSolverHPC() { + delete _mesh; + delete _problem; +} + +void SNSolverHPC::Solve() { + + // --- Preprocessing --- + PrepareVolumeOutput(); + + DrawPreSolverOutput(); + + // Create Backup solution for Runge Kutta + // std::vector solRK0 = _sol; + + auto start = std::chrono::high_resolution_clock::now(); // Start timing + + std::chrono::duration duration; + // Loop over energies (pseudo-time of continuous slowing down approach) + for( unsigned iter = 0; iter < _nIter; iter++ ) { + if( iter == _nIter - 1 ) { // last iteration + _dT = _settings->GetTEnd() - iter * _dT; + // if( _temporalOrder == 2 ) _dT /= 2.0; + } + if( _temporalOrder == 2 ) { + std::vector solRK0( _sol ); + ( _spatialOrder == 2 ) ? FluxOrder2() : FluxOrder1(); + FVMUpdate(); + ( _spatialOrder == 2 ) ? FluxOrder2() : FluxOrder1(); + FVMUpdate(); +#pragma omp parallel for + for( unsigned i = 0; i < _nCells; ++i ) { + _sol[i] = 0.5 * ( solRK0[i] + _sol[i] ); // Solution averaging with HEUN + } + } + else { + ( _spatialOrder == 2 ) ? FluxOrder2() : FluxOrder1(); + FVMUpdate(); + } + + // --- Wall time measurement + duration = std::chrono::high_resolution_clock::now() - start; + _currTime = std::chrono::duration_cast>( duration ).count(); + + // --- Runge Kutta Timestep --- + // if( _temporalOrder == 2 ) RKUpdate( solRK0, _sol ); + // --- Write Output --- + WriteVolumeOutput( iter ); + WriteScalarOutput( iter ); + + // --- Update Scalar Fluxes + + // --- Print Output --- + PrintScreenOutput( iter ); + PrintHistoryOutput( iter ); + PrintVolumeOutput( iter ); + } + + // --- Postprocessing --- + + DrawPostSolverOutput(); +} + +void SNSolverHPC::RKUpdate( std::vector& sol_0, std::vector& sol_rk ) { +#pragma omp parallel for + for( unsigned i = 0; i < _nCells * _nSys; ++i ) { + _sol[i] = 0.5 * ( sol_0[i] + sol_rk[i] ); + } +} + +void SNSolverHPC::PrintVolumeOutput() const { ExportVTK( _settings->GetOutputFile(), _outputFields, _outputFieldNames, _mesh ); } + +void SNSolverHPC::PrintVolumeOutput( int idx_iter ) const { + if( _settings->GetVolumeOutputFrequency() != 0 && idx_iter % (unsigned)_settings->GetVolumeOutputFrequency() == 0 ) { + ExportVTK( _settings->GetOutputFile() + "_" + std::to_string( idx_iter ), _outputFields, _outputFieldNames, _mesh ); // slow + } + if( idx_iter == (int)_nIter - 1 ) { // Last iteration write without suffix. + ExportVTK( _settings->GetOutputFile(), _outputFields, _outputFieldNames, _mesh ); + } +} + +void SNSolverHPC::FluxOrder2() { + + double const eps = 1e-10; + +#pragma omp parallel for + for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) { // Second order information for Flux + unsigned idx_nbr_glob = 0; + if( _cellBoundaryTypes[idx_cell] == BOUNDARY_TYPE::NONE ) { // skip ghost cells + + // #pragma omp simd + for( unsigned idx_sys = 0; idx_sys < _nSys; idx_sys++ ) { + double r; + _limiter[Idx2D( idx_cell, idx_sys, _nSys )] = 1.; // limiter should be zero at boundary + _solDx[Idx3D( idx_cell, idx_sys, 0, _nSys, _nDim )] = 0.; + _solDx[Idx3D( idx_cell, idx_sys, 1, _nSys, _nDim )] = 0.; + + double minSol = _sol[Idx2D( idx_cell, idx_sys, _nSys )]; + double maxSol = minSol; + double solInterfaceAvg = 0.0; + double gaussPoint = 0; + for( unsigned idx_nbr = 0; idx_nbr < _nNbr; ++idx_nbr ) { // Compute slopes and mininum and maximum + idx_nbr_glob = _neighbors[Idx2D( idx_cell, idx_nbr, _nNbr )]; + + // Slopes + solInterfaceAvg = 0.5 * ( _sol[Idx2D( idx_cell, idx_sys, _nSys )] + _sol[Idx2D( idx_nbr_glob, idx_sys, _nSys )] ); + _solDx[Idx3D( idx_cell, idx_sys, 0, _nSys, _nDim )] += solInterfaceAvg * _normals[Idx3D( idx_cell, idx_nbr, 0, _nNbr, _nDim )]; + _solDx[Idx3D( idx_cell, idx_sys, 1, _nSys, _nDim )] += solInterfaceAvg * _normals[Idx3D( idx_cell, idx_nbr, 1, _nNbr, _nDim )]; + + // Compute ptswise max and minimum solultion values of current and neighbor cells + maxSol = std::max( _sol[Idx2D( idx_nbr_glob, idx_sys, _nSys )], maxSol ); + minSol = std::min( _sol[Idx2D( idx_nbr_glob, idx_sys, _nSys )], minSol ); + } + _solDx[Idx3D( idx_cell, idx_sys, 0, _nSys, _nDim )] /= _areas[idx_cell]; + _solDx[Idx3D( idx_cell, idx_sys, 1, _nSys, _nDim )] /= _areas[idx_cell]; + + for( unsigned idx_nbr = 0; idx_nbr < _nNbr; idx_nbr++ ) { // Compute limiter, see https://arxiv.org/pdf/1710.07187.pdf + + // Compute test value at cell vertex, called gaussPt + gaussPoint = + _solDx[Idx3D( idx_cell, idx_sys, 0, _nSys, _nDim )] * _relativeCellVertices[Idx3D( idx_cell, idx_nbr, 0, _nNbr, _nDim )] + + _solDx[Idx3D( idx_cell, idx_sys, 1, _nSys, _nDim )] * _relativeCellVertices[Idx3D( idx_cell, idx_nbr, 1, _nNbr, _nDim )]; + + // BARTH-JESPERSEN LIMITER + r = ( gaussPoint > 0 ) ? std::min( ( maxSol - _sol[Idx2D( idx_cell, idx_sys, _nSys )] ) / ( gaussPoint + eps ), 1.0 ) + : std::min( ( minSol - _sol[Idx2D( idx_cell, idx_sys, _nSys )] ) / ( gaussPoint - eps ), 1.0 ); + + r = ( std::abs( gaussPoint ) < eps ) ? 1 : r; + + _limiter[Idx2D( idx_cell, idx_sys, _nSys )] = std::min( r, _limiter[Idx2D( idx_cell, idx_sys, _nSys )] ); + + // VENKATAKRISHNAN LIMITER + // double delta1Max = maxSol - _sol[Idx2D( idx_cell, idx_sys, _nSys )]; + // double delta1Min = minSol - _sol[Idx2D( idx_cell, idx_sys, _nSys )]; + // + // r = ( gaussPoint > 0 ) ? ( ( delta1Max * delta1Max + _areas[idx_cell] ) * gaussPoint + 2 * gaussPoint * gaussPoint * delta1Max + // ) / + // ( delta1Max * delta1Max + 2 * gaussPoint * gaussPoint + delta1Max * gaussPoint + _areas[idx_cell] ) + // / ( gaussPoint + eps ) + // : ( ( delta1Min * delta1Min + _areas[idx_cell] ) * gaussPoint + 2 * gaussPoint * gaussPoint * delta1Min ) + // / + // ( delta1Min * delta1Min + 2 * gaussPoint * gaussPoint + delta1Min * gaussPoint + _areas[idx_cell] ) + // / ( gaussPoint - eps ); + // + // r = ( std::abs( gaussPoint ) < eps ) ? 1 : r; + // + //_limiter[Idx2D( idx_cell, idx_sys, _nSys )] = std::min( r, _limiter[Idx2D( idx_cell, idx_sys, _nSys )] ); + } + } + } + else { +#pragma omp simd + for( unsigned idx_sys = 0; idx_sys < _nSys; idx_sys++ ) { + _limiter[Idx2D( idx_cell, idx_sys, _nSys )] = 0.; // limiter should be zero at boundary + _solDx[Idx3D( idx_cell, idx_sys, 0, _nSys, _nDim )] = 0.; + _solDx[Idx3D( idx_cell, idx_sys, 1, _nSys, _nDim )] = 0.; + } + } + } + +#pragma omp parallel for + for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) { + +#pragma omp simd + for( unsigned idx_sys = 0; idx_sys < _nSys; idx_sys++ ) { + _flux[Idx2D( idx_cell, idx_sys, _nSys )] = 0.; + } + + // Fluxes + for( unsigned idx_nbr = 0; idx_nbr < _nNbr; ++idx_nbr ) { + if( _cellBoundaryTypes[idx_cell] == BOUNDARY_TYPE::NEUMANN && _neighbors[Idx2D( idx_cell, idx_nbr, _nNbr )] == _nCells ) { +#pragma omp simd + for( unsigned idx_sys = 0; idx_sys < _nSys; idx_sys++ ) { + double localInner = _quadPts[Idx2D( idx_sys, 0, _nDim )] * _normals[Idx3D( idx_cell, idx_nbr, 0, _nNbr, _nDim )] + + _quadPts[Idx2D( idx_sys, 1, _nDim )] * _normals[Idx3D( idx_cell, idx_nbr, 1, _nNbr, _nDim )]; + if( localInner > 0 ) { + _flux[Idx2D( idx_cell, idx_sys, _nSys )] += localInner * _sol[Idx2D( idx_cell, idx_sys, _nSys )]; + } + else { + double ghostCellValue = ( _ghostCellsReflectingY[idx_cell] ) + ? 0.0 * _sol[Idx2D( idx_cell, _quadratureYReflection[idx_sys], _nSys )] // Relecting boundary + : _ghostCells[idx_cell][idx_sys]; // fixed boundary + _flux[Idx2D( idx_cell, idx_sys, _nSys )] += localInner * ghostCellValue; + } + } + } + else { +// Second order +#pragma omp simd + for( unsigned idx_sys = 0; idx_sys < _nSys; idx_sys++ ) { + // store flux contribution on psiNew_sigmaS to save memory + unsigned idx_nbr_glob = _neighbors[Idx2D( idx_cell, idx_nbr, _nNbr )]; // global idx of neighbor cell + double localInner = _quadPts[Idx2D( idx_sys, 0, _nDim )] * _normals[Idx3D( idx_cell, idx_nbr, 0, _nNbr, _nDim )] + + _quadPts[Idx2D( idx_sys, 1, _nDim )] * _normals[Idx3D( idx_cell, idx_nbr, 1, _nNbr, _nDim )]; + + _flux[Idx2D( idx_cell, idx_sys, _nSys )] += + ( localInner > 0 ) ? localInner * ( _sol[Idx2D( idx_cell, idx_sys, _nSys )] + + _limiter[Idx2D( idx_cell, idx_sys, _nSys )] * + ( _solDx[Idx3D( idx_cell, idx_sys, 0, _nSys, _nDim )] * + _relativeInterfaceMidPt[Idx3D( idx_cell, idx_nbr, 0, _nNbr, _nDim )] + + _solDx[Idx3D( idx_cell, idx_sys, 1, _nSys, _nDim )] * + _relativeInterfaceMidPt[Idx3D( idx_cell, idx_nbr, 1, _nNbr, _nDim )] ) ) + : localInner * ( _sol[Idx2D( idx_nbr_glob, idx_sys, _nSys )] + + _limiter[Idx2D( idx_nbr_glob, idx_sys, _nSys )] * + ( _solDx[Idx3D( idx_nbr_glob, idx_sys, 0, _nSys, _nDim )] * + ( _interfaceMidPoints[Idx3D( idx_cell, idx_nbr, 0, _nNbr, _nDim )] - + _cellMidPoints[Idx2D( idx_nbr_glob, 0, _nDim )] ) + + _solDx[Idx3D( idx_nbr_glob, idx_sys, 1, _nSys, _nDim )] * + ( _interfaceMidPoints[Idx3D( idx_cell, idx_nbr, 1, _nNbr, _nDim )] - + _cellMidPoints[Idx2D( idx_nbr_glob, 1, _nDim )] ) ) ); + } + } + } + } +} + +void SNSolverHPC::FluxOrder1() { + +#pragma omp parallel for + for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) { + +#pragma omp simd + for( unsigned idx_sys = 0; idx_sys < _nSys; idx_sys++ ) { + _flux[Idx2D( idx_cell, idx_sys, _nSys )] = 0.0; // Reset temporary variable + } + + // Fluxes + for( unsigned idx_nbr = 0; idx_nbr < _nNbr; ++idx_nbr ) { + if( _cellBoundaryTypes[idx_cell] == BOUNDARY_TYPE::NEUMANN && _neighbors[Idx2D( idx_cell, idx_nbr, _nNbr )] == _nCells ) { +#pragma omp simd + for( unsigned idx_sys = 0; idx_sys < _nSys; idx_sys++ ) { + double localInner = _quadPts[Idx2D( idx_sys, 0, _nDim )] * _normals[Idx3D( idx_cell, idx_nbr, 0, _nNbr, _nDim )] + + _quadPts[Idx2D( idx_sys, 1, _nDim )] * _normals[Idx3D( idx_cell, idx_nbr, 1, _nNbr, _nDim )]; + if( localInner > 0 ) { + _flux[Idx2D( idx_cell, idx_sys, _nSys )] += localInner * _sol[Idx2D( idx_cell, idx_sys, _nSys )]; + } + else { + double ghostCellValue = ( _ghostCellsReflectingY[idx_cell] ) + ? _sol[Idx2D( idx_cell, _quadratureYReflection[idx_sys], _nSys )] // Relecting boundary + : _ghostCells[idx_cell][idx_sys]; // fixed boundary + _flux[Idx2D( idx_cell, idx_sys, _nSys )] += localInner * ghostCellValue; + } + } + } + else { + unsigned nbr_glob = _neighbors[Idx2D( idx_cell, idx_nbr, _nNbr )]; // global idx of neighbor cell +#pragma omp simd + for( unsigned idx_sys = 0; idx_sys < _nSys; idx_sys++ ) { + + double localInner = _quadPts[Idx2D( idx_sys, 0, _nDim )] * _normals[Idx3D( idx_cell, idx_nbr, 0, _nNbr, _nDim )] + + _quadPts[Idx2D( idx_sys, 1, _nDim )] * _normals[Idx3D( idx_cell, idx_nbr, 1, _nNbr, _nDim )]; + + _flux[Idx2D( idx_cell, idx_sys, _nSys )] += ( localInner > 0 ) ? localInner * _sol[Idx2D( idx_cell, idx_sys, _nSys )] + : localInner * _sol[Idx2D( nbr_glob, idx_sys, _nSys )]; + } + } + } + } +} + +void SNSolverHPC::FVMUpdate() { + _mass = 0.0; + _rmsFlux = 0.0; + _curAbsorptionLattice = 0.0; + _curScalarOutflow = 0.0; + _curAbsorptionHohlraumCenter = 0.0; // Green and blue areas of symmetric hohlraum + _curAbsorptionHohlraumVertical = 0.0; // Red areas of symmetric hohlraum + _curAbsorptionHohlraumHorizontal = 0.0; // Black areas of symmetric hohlraum + _varAbsorptionHohlraumGreen = 0.0; + double a_g = 0.0; + +#pragma omp parallel for reduction( + : _mass, \ + _rmsFlux, \ + _curAbsorptionLattice, \ + _curScalarOutflow, \ + _curAbsorptionHohlraumCenter, \ + _curAbsorptionHohlraumVertical, \ + _curAbsorptionHohlraumHorizontal, \ + a_g ) reduction( max : _curMaxOrdinateOutflow, _curMaxAbsorptionLattice ) + for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) { + +#pragma omp simd + for( unsigned idx_sys = 0; idx_sys < _nSys; idx_sys++ ) { + // Upate + _sol[Idx2D( idx_cell, idx_sys, _nSys )] = + ( 1 - _dT * _sigmaT[idx_cell] ) * _sol[Idx2D( idx_cell, idx_sys, _nSys )] - + _dT / _areas[idx_cell] * _flux[Idx2D( idx_cell, idx_sys, _nSys )] + + _dT * ( _sigmaS[idx_cell] * _scalarFlux[idx_cell] / ( 4 * M_PI ) + _source[Idx2D( idx_cell, idx_sys, _nSys )] ); + } + + // --- Iter Postprocessing --- + + double localScalarFlux = 0; + +#pragma omp simd reduction( + : localScalarFlux ) + for( unsigned idx_sys = 0; idx_sys < _nSys; ++idx_sys ) { + _sol[Idx2D( idx_cell, idx_sys, _nSys )] = std::max( _sol[Idx2D( idx_cell, idx_sys, _nSys )], 0.0 ); + localScalarFlux += _sol[Idx2D( idx_cell, idx_sys, _nSys )] * _quadWeights[idx_sys]; + } + + _mass += localScalarFlux * _areas[idx_cell]; + _rmsFlux += ( localScalarFlux - _scalarFlux[idx_cell] ) * ( localScalarFlux - _scalarFlux[idx_cell] ); + _scalarFlux[idx_cell] = localScalarFlux; // set flux + + if( _settings->GetProblemName() == PROBLEM_Lattice || _settings->GetProblemName() == PROBLEM_HalfLattice ) { + if( IsAbsorptionLattice( _cellMidPoints[Idx2D( idx_cell, 0, _nDim )], _cellMidPoints[Idx2D( idx_cell, 1, _nDim )] ) ) { + double sigmaAPsi = _scalarFlux[idx_cell] * ( _sigmaT[idx_cell] - _sigmaS[idx_cell] ); + _curAbsorptionLattice += sigmaAPsi * _areas[idx_cell]; + _curMaxAbsorptionLattice = ( _curMaxAbsorptionLattice < sigmaAPsi ) ? sigmaAPsi : _curMaxAbsorptionLattice; + } + } + + if( _settings->GetProblemName() == PROBLEM_SymmetricHohlraum || _settings->GetProblemName() == PROBLEM_QuarterHohlraum ) { + + double x = _cellMidPoints[Idx2D( idx_cell, 0, _nDim )]; + double y = _cellMidPoints[Idx2D( idx_cell, 1, _nDim )]; + + if( x > -0.2 && x < 0.2 && y > -0.35 && y < 0.35 ) { + _curAbsorptionHohlraumCenter += _scalarFlux[idx_cell] * ( _sigmaT[idx_cell] - _sigmaS[idx_cell] ) * _areas[idx_cell]; + } + if( ( x < -0.6 && y > -0.4 && y < 0.4 ) || ( x > 0.6 && y > -0.4 && y < 0.4 ) ) { + _curAbsorptionHohlraumVertical += _scalarFlux[idx_cell] * ( _sigmaT[idx_cell] - _sigmaS[idx_cell] ) * _areas[idx_cell]; + } + if( y > 0.6 || y < -0.6 ) { + _curAbsorptionHohlraumHorizontal += _scalarFlux[idx_cell] * ( _sigmaT[idx_cell] - _sigmaS[idx_cell] ) * _areas[idx_cell]; + } + + // Variation in absorption of center (part 1) + bool green1 = x > -0.2 && x < -0.15 && y > -0.35 && y < 0.35; // green area 1 (lower boundary) + bool green2 = x > 0.15 && x < 0.2 && y > -0.35 && y < 0.35; // green area 2 (upper boundary) + bool green3 = x > -0.2 && x < 0.2 && y > -0.4 && y < -0.35; // green area 3 (left boundary) + bool green4 = x > -0.2 && x < 0.2 && y > 0.35 && y < 0.4; // green area 4 (right boundary) + + if( green1 || green2 || green3 || green4 ) { + a_g += ( _sigmaT[idx_cell] - _sigmaS[idx_cell] ) * _scalarFlux[idx_cell] * _areas[idx_cell]; + } + } + + // Outflow + if( _cellBoundaryTypes[idx_cell] == BOUNDARY_TYPE::NEUMANN && !_ghostCellsReflectingY[idx_cell] ) { + // Iterate over face cell faces + double currOrdinatewiseOutflow = 0.0; + + for( unsigned idx_nbr = 0; idx_nbr < _nNbr; ++idx_nbr ) { + // Find face that points outward + if( _neighbors[Idx2D( idx_cell, idx_nbr, _nNbr )] == _nCells ) { +#pragma omp simd reduction( + : _curScalarOutflow ) reduction( max : _curMaxOrdinateOutflow ) + for( unsigned idx_sys = 0; idx_sys < _nSys; ++idx_sys ) { + double localInner = _quadPts[Idx2D( idx_sys, 0, _nDim )] * _normals[Idx3D( idx_cell, idx_nbr, 0, _nNbr, _nDim )] + + _quadPts[Idx2D( idx_sys, 1, _nDim )] * _normals[Idx3D( idx_cell, idx_nbr, 1, _nNbr, _nDim )]; + // Find outward facing transport directions + + if( localInner > 0.0 ) { + _curScalarOutflow += localInner * _sol[Idx2D( idx_cell, idx_sys, _nSys )] * _quadWeights[idx_sys]; // Integrate flux + + currOrdinatewiseOutflow = + _sol[Idx2D( idx_cell, idx_sys, _nSys )] * localInner / + sqrt( ( + _normals[Idx3D( idx_cell, idx_nbr, 0, _nNbr, _nDim )] * _normals[Idx3D( idx_cell, idx_nbr, 0, _nNbr, _nDim )] + + _normals[Idx3D( idx_cell, idx_nbr, 1, _nNbr, _nDim )] * _normals[Idx3D( idx_cell, idx_nbr, 1, _nNbr, _nDim )] ) ); + + _curMaxOrdinateOutflow = + ( currOrdinatewiseOutflow > _curMaxOrdinateOutflow ) ? currOrdinatewiseOutflow : _curMaxOrdinateOutflow; + } + } + } + } + } + } + // Variation absorption (part II) + if( _settings->GetProblemName() == PROBLEM_SymmetricHohlraum || _settings->GetProblemName() == PROBLEM_QuarterHohlraum ) { +#pragma omp parallel for reduction( + : _varAbsorptionHohlraumGreen ) + for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) { + double x = _cellMidPoints[Idx2D( idx_cell, 0, _nDim )]; + double y = _cellMidPoints[Idx2D( idx_cell, 1, _nDim )]; + bool green1, green2, green3, green4; + + green1 = x > -0.2 && x < -0.15 && y > -0.35 && y < 0.35; // green area 1 (lower boundary) + green2 = x > 0.15 && x < 0.2 && y > -0.35 && y < 0.35; // green area 2 (upper boundary) + green3 = x > -0.2 && x < 0.2 && y > -0.4 && y < -0.35; // green area 3 (left boundary) + green4 = x > -0.2 && x < 0.2 && y > 0.35 && y < 0.4; // green area 4 (right boundary) + + if( green1 || green2 || green3 || green4 ) { + _varAbsorptionHohlraumGreen += ( a_g - _scalarFlux[idx_cell] * ( _sigmaT[idx_cell] - _sigmaS[idx_cell] ) ) * + ( a_g - _scalarFlux[idx_cell] * ( _sigmaT[idx_cell] - _sigmaS[idx_cell] ) ) * _areas[idx_cell]; + } + } +// Probes value moments +#pragma omp parallel for + for( unsigned idx_probe = 0; idx_probe < 4; idx_probe++ ) { // Loop over probing cells + _probingMoments[Idx2D( idx_probe, 0, 3 )] = _scalarFlux[_probingCellsHohlraum[idx_probe]]; + _probingMoments[Idx2D( idx_probe, 1, 3 )] = 0.0; + _probingMoments[Idx2D( idx_probe, 2, 3 )] = 0.0; + + // for( unsigned idx_sys = 0; idx_sys < _nSys; idx_sys++ ) { + // _probingMoments[Idx2D( idx_probe, 1, 3 )] += + // _quadPts[Idx2D( idx_sys, 0, _nDim )] * _sol[Idx2D( _probingCellsHohlraum[idx_probe], idx_sys, _nSys )] * _quadWeights[idx_sys]; + // _probingMoments[Idx2D( idx_probe, 2, 3 )] += + // _quadPts[Idx2D( idx_sys, 1, _nDim )] * _sol[Idx2D( _probingCellsHohlraum[idx_probe], idx_sys, _nSys )] * _quadWeights[idx_sys]; + // } + } + + // probe values green + ComputeQOIsGreenProbingLine(); + } + + _rmsFlux = sqrt( _rmsFlux ); + _totalScalarOutflow += _curScalarOutflow * _dT; + _totalAbsorptionLattice += _curAbsorptionLattice * _dT; + + _totalAbsorptionHohlraumCenter += _curAbsorptionHohlraumCenter * _dT; + _totalAbsorptionHohlraumVertical += _curAbsorptionHohlraumVertical * _dT; + _totalAbsorptionHohlraumHorizontal += _curAbsorptionHohlraumHorizontal * _dT; +} + +bool SNSolverHPC::IsAbsorptionLattice( double x, double y ) const { + // Check whether pos is inside absorbing squares + double xy_corrector = -3.5; + std::vector lbounds{ 1 + xy_corrector, 2 + xy_corrector, 3 + xy_corrector, 4 + xy_corrector, 5 + xy_corrector }; + std::vector ubounds{ 2 + xy_corrector, 3 + xy_corrector, 4 + xy_corrector, 5 + xy_corrector, 6 + xy_corrector }; + for( unsigned k = 0; k < lbounds.size(); ++k ) { + for( unsigned l = 0; l < lbounds.size(); ++l ) { + if( ( l + k ) % 2 == 1 || ( k == 2 && l == 2 ) || ( k == 2 && l == 4 ) ) continue; + if( x >= lbounds[k] && x <= ubounds[k] && y >= lbounds[l] && y <= ubounds[l] ) { + return true; + } + } + } + return false; +} + +// --- Helper --- +double SNSolverHPC::ComputeTimeStep( double cfl ) const { + // for pseudo 1D, set timestep to dx + double dx, dy; + switch( _settings->GetProblemName() ) { + case PROBLEM_Checkerboard1D: + dx = 7.0 / (double)_nCells; + dy = 0.3; + return cfl * ( dx * dy ) / ( dx + dy ); + break; + case PROBLEM_Linesource1D: // Fallthrough + case PROBLEM_Meltingcube1D: // Fallthrough + case PROBLEM_Aircavity1D: + dx = 3.0 / (double)_nCells; + dy = 0.3; + return cfl * ( dx * dy ) / ( dx + dy ); + break; + default: break; // 2d as normal + } + // 2D case + double charSize = __DBL_MAX__; // minimum char size of all mesh cells in the mesh + for( unsigned j = 0; j < _nCells; j++ ) { + double currCharSize = sqrt( _areas[j] ); + if( currCharSize < charSize ) { + charSize = currCharSize; + } + } + auto log = spdlog::get( "event" ); + std::string line = "| Smallest characteristic length of a grid cell in this mesh: " + std::to_string( charSize ); + log->info( line ); + line = "| Corresponding maximal time-step: " + std::to_string( cfl * charSize ); + log->info( line ); + return cfl * charSize; +} + +// --- IO ---- +void SNSolverHPC::PrepareScreenOutput() { + unsigned nFields = (unsigned)_settings->GetNScreenOutput(); + + _screenOutputFieldNames.resize( nFields ); + _screenOutputFields.resize( nFields ); + + // Prepare all output Fields ==> Specified in option SCREEN_OUTPUT + for( unsigned idx_field = 0; idx_field < nFields; idx_field++ ) { + // Prepare all Output Fields per group + + // Different procedure, depending on the Group... + switch( _settings->GetScreenOutput()[idx_field] ) { + case MASS: _screenOutputFieldNames[idx_field] = "Mass"; break; + case ITER: _screenOutputFieldNames[idx_field] = "Iter"; break; + case WALL_TIME: _screenOutputFieldNames[idx_field] = "Wall time [s]"; break; + case RMS_FLUX: _screenOutputFieldNames[idx_field] = "RMS flux"; break; + case VTK_OUTPUT: _screenOutputFieldNames[idx_field] = "VTK out"; break; + case CSV_OUTPUT: _screenOutputFieldNames[idx_field] = "CSV out"; break; + case CUR_OUTFLOW: _screenOutputFieldNames[idx_field] = "Cur. outflow"; break; + case TOTAL_OUTFLOW: _screenOutputFieldNames[idx_field] = "Tot. outflow"; break; + case MAX_OUTFLOW: _screenOutputFieldNames[idx_field] = "Max outflow"; break; + case CUR_PARTICLE_ABSORPTION: _screenOutputFieldNames[idx_field] = "Cur. absorption"; break; + case TOTAL_PARTICLE_ABSORPTION: _screenOutputFieldNames[idx_field] = "Tot. absorption"; break; + case MAX_PARTICLE_ABSORPTION: _screenOutputFieldNames[idx_field] = "Max absorption"; break; + case TOTAL_PARTICLE_ABSORPTION_CENTER: _screenOutputFieldNames[idx_field] = "Tot. abs. center"; break; + case TOTAL_PARTICLE_ABSORPTION_VERTICAL: _screenOutputFieldNames[idx_field] = "Tot. abs. vertical wall"; break; + case TOTAL_PARTICLE_ABSORPTION_HORIZONTAL: _screenOutputFieldNames[idx_field] = "Tot. abs. horizontal wall"; break; + case PROBE_MOMENT_TIME_TRACE: + _screenOutputFieldNames[idx_field] = "Probe 1 u_0"; + idx_field++; + _screenOutputFieldNames[idx_field] = "Probe 2 u_0"; + if( _settings->GetProblemName() == PROBLEM_SymmetricHohlraum ) { + idx_field++; + _screenOutputFieldNames[idx_field] = "Probe 3 u_0"; + idx_field++; + _screenOutputFieldNames[idx_field] = "Probe 4 u_0"; + } + break; + case VAR_ABSORPTION_GREEN: _screenOutputFieldNames[idx_field] = "Var. absorption green"; break; + default: ErrorMessages::Error( "Screen output field not defined!", CURRENT_FUNCTION ); break; + } + } +} + +void SNSolverHPC::WriteScalarOutput( unsigned idx_iter ) { + unsigned n_probes = 4; + + unsigned nFields = (unsigned)_settings->GetNScreenOutput(); + const VectorVector probingMoments = _problem->GetCurrentProbeMoment(); + // -- Screen Output + for( unsigned idx_field = 0; idx_field < nFields; idx_field++ ) { + // Prepare all Output Fields per group + // Different procedure, depending on the Group... + switch( _settings->GetScreenOutput()[idx_field] ) { + case MASS: _screenOutputFields[idx_field] = _mass; break; + case ITER: _screenOutputFields[idx_field] = idx_iter; break; + case WALL_TIME: _screenOutputFields[idx_field] = _currTime; break; + case RMS_FLUX: _screenOutputFields[idx_field] = _rmsFlux; break; + case VTK_OUTPUT: + _screenOutputFields[idx_field] = 0; + if( ( _settings->GetVolumeOutputFrequency() != 0 && idx_iter % (unsigned)_settings->GetVolumeOutputFrequency() == 0 ) || + ( idx_iter == _nIter - 1 ) /* need sol at last iteration */ ) { + _screenOutputFields[idx_field] = 1; + } + break; + case CSV_OUTPUT: + _screenOutputFields[idx_field] = 0; + if( ( _settings->GetHistoryOutputFrequency() != 0 && idx_iter % (unsigned)_settings->GetHistoryOutputFrequency() == 0 ) || + ( idx_iter == _nIter - 1 ) /* need sol at last iteration */ ) { + _screenOutputFields[idx_field] = 1; + } + break; + case CUR_OUTFLOW: _screenOutputFields[idx_field] = _curScalarOutflow; break; + case TOTAL_OUTFLOW: _screenOutputFields[idx_field] = _totalScalarOutflow; break; + case MAX_OUTFLOW: _screenOutputFields[idx_field] = _curMaxOrdinateOutflow; break; + case CUR_PARTICLE_ABSORPTION: _screenOutputFields[idx_field] = _curAbsorptionLattice; break; + case TOTAL_PARTICLE_ABSORPTION: _screenOutputFields[idx_field] = _totalAbsorptionLattice; break; + case MAX_PARTICLE_ABSORPTION: _screenOutputFields[idx_field] = _curMaxAbsorptionLattice; break; + case TOTAL_PARTICLE_ABSORPTION_CENTER: _screenOutputFields[idx_field] = _totalAbsorptionHohlraumCenter; break; + case TOTAL_PARTICLE_ABSORPTION_VERTICAL: _screenOutputFields[idx_field] = _totalAbsorptionHohlraumVertical; break; + case TOTAL_PARTICLE_ABSORPTION_HORIZONTAL: _screenOutputFields[idx_field] = _totalAbsorptionHohlraumHorizontal; break; + case PROBE_MOMENT_TIME_TRACE: + if( _settings->GetProblemName() == PROBLEM_SymmetricHohlraum ) n_probes = 4; + if( _settings->GetProblemName() == PROBLEM_QuarterHohlraum ) n_probes = 2; + for( unsigned i = 0; i < n_probes; i++ ) { + _screenOutputFields[idx_field] = _probingMoments[Idx2D( i, 0, 3 )]; + idx_field++; + } + idx_field--; + break; + case VAR_ABSORPTION_GREEN: _screenOutputFields[idx_field] = _varAbsorptionHohlraumGreen; break; + default: ErrorMessages::Error( "Screen output group not defined!", CURRENT_FUNCTION ); break; + } + } + + // --- History output --- + nFields = (unsigned)_settings->GetNHistoryOutput(); + + std::vector screenOutputFields = _settings->GetScreenOutput(); + for( unsigned idx_field = 0; idx_field < nFields; idx_field++ ) { + + // Prepare all Output Fields per group + // Different procedure, depending on the Group... + switch( _settings->GetHistoryOutput()[idx_field] ) { + case MASS: _historyOutputFields[idx_field] = _mass; break; + case ITER: _historyOutputFields[idx_field] = idx_iter; break; + case WALL_TIME: _historyOutputFields[idx_field] = _currTime; break; + case RMS_FLUX: _historyOutputFields[idx_field] = _rmsFlux; break; + case VTK_OUTPUT: + _historyOutputFields[idx_field] = 0; + if( ( _settings->GetVolumeOutputFrequency() != 0 && idx_iter % (unsigned)_settings->GetVolumeOutputFrequency() == 0 ) || + ( idx_iter == _nIter - 1 ) /* need sol at last iteration */ ) { + _historyOutputFields[idx_field] = 1; + } + break; + + case CSV_OUTPUT: + _historyOutputFields[idx_field] = 0; + if( ( _settings->GetHistoryOutputFrequency() != 0 && idx_iter % (unsigned)_settings->GetHistoryOutputFrequency() == 0 ) || + ( idx_iter == _nIter - 1 ) /* need sol at last iteration */ ) { + _historyOutputFields[idx_field] = 1; + } + break; + case CUR_OUTFLOW: _historyOutputFields[idx_field] = _curScalarOutflow; break; + case TOTAL_OUTFLOW: _historyOutputFields[idx_field] = _totalScalarOutflow; break; + case MAX_OUTFLOW: _historyOutputFields[idx_field] = _curMaxOrdinateOutflow; break; + case CUR_PARTICLE_ABSORPTION: _historyOutputFields[idx_field] = _curAbsorptionLattice; break; + case TOTAL_PARTICLE_ABSORPTION: _historyOutputFields[idx_field] = _totalAbsorptionLattice; break; + case MAX_PARTICLE_ABSORPTION: _historyOutputFields[idx_field] = _curMaxAbsorptionLattice; break; + case TOTAL_PARTICLE_ABSORPTION_CENTER: _historyOutputFields[idx_field] = _totalAbsorptionHohlraumCenter; break; + case TOTAL_PARTICLE_ABSORPTION_VERTICAL: _historyOutputFields[idx_field] = _totalAbsorptionHohlraumVertical; break; + case TOTAL_PARTICLE_ABSORPTION_HORIZONTAL: _historyOutputFields[idx_field] = _totalAbsorptionHohlraumHorizontal; break; + case PROBE_MOMENT_TIME_TRACE: + if( _settings->GetProblemName() == PROBLEM_SymmetricHohlraum ) n_probes = 4; + if( _settings->GetProblemName() == PROBLEM_QuarterHohlraum ) n_probes = 2; + for( unsigned i = 0; i < n_probes; i++ ) { + for( unsigned j = 0; j < 3; j++ ) { + _screenOutputFields[idx_field] = _probingMoments[Idx2D( i, j, 3 )]; + idx_field++; + } + } + idx_field--; + break; + case VAR_ABSORPTION_GREEN: _historyOutputFields[idx_field] = _varAbsorptionHohlraumGreen; break; + case VAR_ABSORPTION_GREEN_LINE: + for( unsigned i = 0; i < _settings->GetNumProbingCellsLineHohlraum(); i++ ) { + _historyOutputFieldNames[idx_field] = _varAbsorptionValsIntegrated[i]; + idx_field++; + } + idx_field--; + break; + default: ErrorMessages::Error( "History output group not defined!", CURRENT_FUNCTION ); break; + } + } +} + +void SNSolverHPC::PrintScreenOutput( unsigned idx_iter ) { + auto log = spdlog::get( "event" ); + + unsigned strLen = 15; // max width of one column + char paddingChar = ' '; + + // assemble the line to print + std::string lineToPrint = "| "; + std::string tmp; + for( unsigned idx_field = 0; idx_field < _settings->GetNScreenOutput(); idx_field++ ) { + tmp = std::to_string( _screenOutputFields[idx_field] ); + + // Format outputs correctly + std::vector integerFields = { ITER }; + std::vector scientificFields = { RMS_FLUX, + MASS, + WALL_TIME, + CUR_OUTFLOW, + TOTAL_OUTFLOW, + MAX_OUTFLOW, + CUR_PARTICLE_ABSORPTION, + TOTAL_PARTICLE_ABSORPTION, + MAX_PARTICLE_ABSORPTION, + TOTAL_PARTICLE_ABSORPTION_CENTER, + TOTAL_PARTICLE_ABSORPTION_VERTICAL, + TOTAL_PARTICLE_ABSORPTION_HORIZONTAL, + PROBE_MOMENT_TIME_TRACE, + VAR_ABSORPTION_GREEN, + VAR_ABSORPTION_GREEN_LINE }; + std::vector booleanFields = { VTK_OUTPUT, CSV_OUTPUT }; + + if( !( integerFields.end() == std::find( integerFields.begin(), integerFields.end(), _settings->GetScreenOutput()[idx_field] ) ) ) { + tmp = std::to_string( (int)_screenOutputFields[idx_field] ); + } + else if( !( booleanFields.end() == std::find( booleanFields.begin(), booleanFields.end(), _settings->GetScreenOutput()[idx_field] ) ) ) { + tmp = "no"; + if( (bool)_screenOutputFields[idx_field] ) tmp = "yes"; + } + else if( !( scientificFields.end() == + std::find( scientificFields.begin(), scientificFields.end(), _settings->GetScreenOutput()[idx_field] ) ) ) { + + std::stringstream ss; + ss << TextProcessingToolbox::DoubleToScientificNotation( _screenOutputFields[idx_field] ); + tmp = ss.str(); + tmp.erase( std::remove( tmp.begin(), tmp.end(), '+' ), tmp.end() ); // removing the '+' sign + } + + if( strLen > tmp.size() ) // Padding + tmp.insert( 0, strLen - tmp.size(), paddingChar ); + else if( strLen < tmp.size() ) // Cutting + tmp.resize( strLen ); + + lineToPrint += tmp + " |"; + } + if( _settings->GetScreenOutputFrequency() != 0 && idx_iter % (unsigned)_settings->GetScreenOutputFrequency() == 0 ) { + log->info( lineToPrint ); + } + else if( idx_iter == _nIter - 1 ) { // Always print last iteration + log->info( lineToPrint ); + } +} + +void SNSolverHPC::PrepareHistoryOutput() { + unsigned n_probes = 4; + + unsigned nFields = (unsigned)_settings->GetNHistoryOutput(); + + _historyOutputFieldNames.resize( nFields ); + _historyOutputFields.resize( nFields ); + + // Prepare all output Fields ==> Specified in option SCREEN_OUTPUT + for( unsigned idx_field = 0; idx_field < nFields; idx_field++ ) { + // Prepare all Output Fields per group + + // Different procedure, depending on the Group... + switch( _settings->GetHistoryOutput()[idx_field] ) { + case MASS: _historyOutputFieldNames[idx_field] = "Mass"; break; + case ITER: _historyOutputFieldNames[idx_field] = "Iter"; break; + case WALL_TIME: _historyOutputFieldNames[idx_field] = "Wall_time_[s]"; break; + case RMS_FLUX: _historyOutputFieldNames[idx_field] = "RMS_flux"; break; + case VTK_OUTPUT: _historyOutputFieldNames[idx_field] = "VTK_out"; break; + case CSV_OUTPUT: _historyOutputFieldNames[idx_field] = "CSV_out"; break; + case CUR_OUTFLOW: _historyOutputFieldNames[idx_field] = "Cur_outflow"; break; + case TOTAL_OUTFLOW: _historyOutputFieldNames[idx_field] = "Total_outflow"; break; + case MAX_OUTFLOW: _historyOutputFieldNames[idx_field] = "Max_outflow"; break; + case CUR_PARTICLE_ABSORPTION: _historyOutputFieldNames[idx_field] = "Cur_absorption"; break; + case TOTAL_PARTICLE_ABSORPTION: _historyOutputFieldNames[idx_field] = "Total_absorption"; break; + case MAX_PARTICLE_ABSORPTION: _historyOutputFieldNames[idx_field] = "Max_absorption"; break; + case TOTAL_PARTICLE_ABSORPTION_CENTER: _historyOutputFieldNames[idx_field] = "Cumulated_absorption_center"; break; + case TOTAL_PARTICLE_ABSORPTION_VERTICAL: _historyOutputFieldNames[idx_field] = "Cumulated_absorption_vertical_wall"; break; + case TOTAL_PARTICLE_ABSORPTION_HORIZONTAL: _historyOutputFieldNames[idx_field] = "Cumulated_absorption_horizontal_wall"; break; + case PROBE_MOMENT_TIME_TRACE: + if( _settings->GetProblemName() == PROBLEM_SymmetricHohlraum ) n_probes = 4; + if( _settings->GetProblemName() == PROBLEM_QuarterHohlraum ) n_probes = 2; + for( unsigned i = 0; i < n_probes; i++ ) { + for( unsigned j = 0; j < 3; j++ ) { + _historyOutputFieldNames[idx_field] = "Probe " + std::to_string( i ) + " u_" + std::to_string( j ); + idx_field++; + } + } + idx_field--; + break; + case VAR_ABSORPTION_GREEN: _historyOutputFieldNames[idx_field] = "Var. absorption green"; break; + case VAR_ABSORPTION_GREEN_LINE: + for( unsigned i = 0; i < _settings->GetNumProbingCellsLineHohlraum(); i++ ) { + _historyOutputFieldNames[idx_field] = "Probe Green Line " + std::to_string( i ); + idx_field++; + } + idx_field--; + break; + default: ErrorMessages::Error( "History output field not defined!", CURRENT_FUNCTION ); break; + } + } +} + +void SNSolverHPC::PrintHistoryOutput( unsigned idx_iter ) { + + auto log = spdlog::get( "tabular" ); + + // assemble the line to print + std::string lineToPrint = ""; + std::string tmp; + for( int idx_field = 0; idx_field < _settings->GetNHistoryOutput() - 1; idx_field++ ) { + if( idx_field == 0 ) { + tmp = std::to_string( _historyOutputFields[idx_field] ); // Iteration count + } + else { + tmp = TextProcessingToolbox::DoubleToScientificNotation( _historyOutputFields[idx_field] ); + } + lineToPrint += tmp + ","; + } + tmp = TextProcessingToolbox::DoubleToScientificNotation( _historyOutputFields[_settings->GetNScreenOutput() - 1] ); + lineToPrint += tmp; // Last element without comma + + if( _settings->GetHistoryOutputFrequency() != 0 && idx_iter % (unsigned)_settings->GetHistoryOutputFrequency() == 0 ) { + log->info( lineToPrint ); + } + else if( idx_iter == _nIter - 1 ) { // Always print last iteration + log->info( lineToPrint ); + } +} + +void SNSolverHPC::DrawPreSolverOutput() { + + // Logger + auto log = spdlog::get( "event" ); + auto logCSV = spdlog::get( "tabular" ); + + std::string hLine = "--"; + + unsigned strLen = 15; // max width of one column + char paddingChar = ' '; + + // Assemble Header for Screen Output + std::string lineToPrint = "| "; + std::string tmpLine = "-----------------"; + for( unsigned idxFields = 0; idxFields < _settings->GetNScreenOutput(); idxFields++ ) { + std::string tmp = _screenOutputFieldNames[idxFields]; + + if( strLen > tmp.size() ) // Padding + tmp.insert( 0, strLen - tmp.size(), paddingChar ); + else if( strLen < tmp.size() ) // Cutting + tmp.resize( strLen ); + + lineToPrint += tmp + " |"; + hLine += tmpLine; + } + log->info( "---------------------------- Solver Starts -----------------------------" ); + log->info( "| The simulation will run for {} iterations.", _nIter ); + log->info( "| The spatial grid contains {} cells.", _nCells ); + if( _settings->GetSolverName() != PN_SOLVER && _settings->GetSolverName() != CSD_PN_SOLVER ) { + log->info( "| The velocity grid contains {} points.", _nq ); + } + log->info( hLine ); + log->info( lineToPrint ); + log->info( hLine ); + + std::string lineToPrintCSV = ""; + for( int idxFields = 0; idxFields < _settings->GetNHistoryOutput() - 1; idxFields++ ) { + std::string tmp = _historyOutputFieldNames[idxFields]; + lineToPrintCSV += tmp + ","; + } + lineToPrintCSV += _historyOutputFieldNames[_settings->GetNHistoryOutput() - 1]; + logCSV->info( lineToPrintCSV ); +} + +void SNSolverHPC::DrawPostSolverOutput() { + + // Logger + auto log = spdlog::get( "event" ); + + std::string hLine = "--"; + + unsigned strLen = 10; // max width of one column + char paddingChar = ' '; + + // Assemble Header for Screen Output + std::string lineToPrint = "| "; + std::string tmpLine = "------------"; + for( unsigned idxFields = 0; idxFields < _settings->GetNScreenOutput(); idxFields++ ) { + std::string tmp = _screenOutputFieldNames[idxFields]; + + if( strLen > tmp.size() ) // Padding + tmp.insert( 0, strLen - tmp.size(), paddingChar ); + else if( strLen < tmp.size() ) // Cutting + tmp.resize( strLen ); + + lineToPrint += tmp + " |"; + hLine += tmpLine; + } + log->info( hLine ); +#ifndef BUILD_TESTING + log->info( "| The volume output files have been stored at " + _settings->GetOutputFile() ); + log->info( "| The log files have been stored at " + _settings->GetLogDir() + _settings->GetLogFile() ); +#endif + log->info( "--------------------------- Solver Finished ----------------------------" ); +} + +unsigned SNSolverHPC::Idx2D( unsigned idx1, unsigned idx2, unsigned len2 ) { return idx1 * len2 + idx2; } + +unsigned SNSolverHPC::Idx3D( unsigned idx1, unsigned idx2, unsigned idx3, unsigned len2, unsigned len3 ) { + return ( idx1 * len2 + idx2 ) * len3 + idx3; +} + +void SNSolverHPC::WriteVolumeOutput( unsigned idx_iter ) { + unsigned nGroups = (unsigned)_settings->GetNVolumeOutput(); + if( ( _settings->GetVolumeOutputFrequency() != 0 && idx_iter % (unsigned)_settings->GetVolumeOutputFrequency() == 0 ) || + ( idx_iter == _nIter - 1 ) /* need sol at last iteration */ ) { + for( unsigned idx_group = 0; idx_group < nGroups; idx_group++ ) { + switch( _settings->GetVolumeOutput()[idx_group] ) { + case MINIMAL: + for( unsigned idx_cell = 0; idx_cell < _nCells; ++idx_cell ) { + _outputFields[idx_group][0][idx_cell] = _scalarFlux[idx_cell]; + } + break; + + default: ErrorMessages::Error( "Volume Output Group not defined for HPC SN Solver!", CURRENT_FUNCTION ); break; + } + } + } +} + +void SNSolverHPC::PrepareVolumeOutput() { + unsigned nGroups = (unsigned)_settings->GetNVolumeOutput(); + + _outputFieldNames.resize( nGroups ); + _outputFields.resize( nGroups ); + + // Prepare all OutputGroups ==> Specified in option VOLUME_OUTPUT + for( unsigned idx_group = 0; idx_group < nGroups; idx_group++ ) { + // Prepare all Output Fields per group + + // Different procedure, depending on the Group... + switch( _settings->GetVolumeOutput()[idx_group] ) { + case MINIMAL: + // Currently only one entry ==> rad flux + _outputFields[idx_group].resize( 1 ); + _outputFieldNames[idx_group].resize( 1 ); + + _outputFields[idx_group][0].resize( _nCells ); + _outputFieldNames[idx_group][0] = "scalar flux"; + break; + + default: ErrorMessages::Error( "Volume Output Group not defined for HPC SN Solver!", CURRENT_FUNCTION ); break; + } + } +} + +void SNSolverHPC::SetGhostCells() { + if( _settings->GetProblemName() == PROBLEM_Lattice ) { + std::vector void_flow( _nSys, 0.0 ); + + for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) { + if( _cellBoundaryTypes[idx_cell] == BOUNDARY_TYPE::NEUMANN || _cellBoundaryTypes[idx_cell] == BOUNDARY_TYPE::DIRICHLET ) { + _ghostCells[idx_cell] = void_flow; + _ghostCellsReflectingY[idx_cell] = false; + } + } + } + else if( _settings->GetProblemName() == PROBLEM_HalfLattice ) { // HALF LATTICE NOT WORKING + + double tol = 1e-12; // For distance to boundary + + if( _settings->GetQuadName() != QUAD_GaussLegendreTensorized2D ) { + ErrorMessages::Error( "This simplified test case only works with symmetric quadrature orders. Use QUAD_GAUSS_LEGENDRE_TENSORIZED_2D", + CURRENT_FUNCTION ); + } + + { // Create the symmetry maps for the quadratures + for( unsigned idx_q = 0; idx_q < _nSys; idx_q++ ) { + for( unsigned idx_q2 = 0; idx_q2 < _nSys; idx_q2++ ) { + if( abs( _quadPts[Idx2D( idx_q, 0, _nDim )] + _quadPts[Idx2D( idx_q2, 0, _nDim )] ) + + abs( _quadPts[Idx2D( idx_q, 1, _nDim )] - _quadPts[Idx2D( idx_q2, 1, _nDim )] ) < + tol ) { + _quadratureYReflection[idx_q] = idx_q2; + break; + } + } + } + } + + if( _quadratureYReflection.size() != _nSys ) { + ErrorMessages::Error( "Problem with Y symmetry of quadrature of this mesh", CURRENT_FUNCTION ); + } + + auto nodes = _mesh->GetNodes(); + auto cellNodes = _mesh->GetCells(); + + for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) { + if( _cellBoundaryTypes[idx_cell] != BOUNDARY_TYPE::NONE ) { + + _ghostCellsReflectingY[idx_cell] = false; + for( unsigned idx_node = 0; idx_node < _nNbr; idx_node++ ) { // Check if corner node is in this cell + if( abs( nodes[cellNodes[idx_cell][idx_node]][1] ) < -3.5 + tol || abs( nodes[cellNodes[idx_cell][idx_node]][1] ) > 3.5 - tol ) { + // upper and lower boundary + _ghostCells.insert( { idx_cell, std::vector( _nSys, 1e-15 ) } ); + break; + } + else if( abs( nodes[cellNodes[idx_cell][idx_node]][0] ) < tol ) { // close to 0 => left boundary + _ghostCellsReflectingY[idx_cell] = true; + // _ghostCells.insert( { idx_cell, std::vector( _nSys ) } ); + break; + } + else { // right boundary + _ghostCells.insert( { idx_cell, std::vector( _nSys, 1e-15 ) } ); + break; + } + } + } + } + } + + else if( _settings->GetProblemName() == PROBLEM_SymmetricHohlraum ) { + std::vector left_inflow( _nSys, 0.0 ); + std::vector right_inflow( _nSys, 0.0 ); + std::vector vertical_flow( _nSys, 0.0 ); + + for( unsigned idx_q = 0; idx_q < _nSys; idx_q++ ) { + if( _quadPts[Idx2D( idx_q, 0, _nDim )] > 0.0 ) left_inflow[idx_q] = 1.0; + if( _quadPts[Idx2D( idx_q, 1, _nDim )] < 0.0 ) right_inflow[idx_q] = 1.0; + } + + for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) { + + _ghostCellsReflectingY[idx_cell] = false; + + if( _cellBoundaryTypes[idx_cell] == BOUNDARY_TYPE::NEUMANN || _cellBoundaryTypes[idx_cell] == BOUNDARY_TYPE::DIRICHLET ) { + if( _cellMidPoints[Idx2D( idx_cell, 1, _nDim )] < -0.6 ) + _ghostCells[idx_cell] = vertical_flow; + else if( _cellMidPoints[Idx2D( idx_cell, 1, _nDim )] > 0.6 ) + _ghostCells[idx_cell] = vertical_flow; + else if( _cellMidPoints[Idx2D( idx_cell, 0, _nDim )] < -0.6 ) + _ghostCells[idx_cell] = left_inflow; + else if( _cellMidPoints[Idx2D( idx_cell, 0, _nDim )] > 0.6 ) + _ghostCells[idx_cell] = right_inflow; + } + } + } +} + +void SNSolverHPC::SetProbingCellsLineGreen() { + std::cout << "SetProbingCellsLineGreen" << std::endl; + + double verticalLineWidth = std::abs( _cornerUpperLeftGreen[1] - _cornerLowerLeftGreen[1] ); + double horizontalLineWidth = std::abs( _cornerUpperLeftGreen[0] - _cornerUpperRightGreen[0] ); + + // double dx = 2 * ( horizontalLineWidth + verticalLineWidth ) / ( (double)_nProbingCellsLineGreen ); + + unsigned nHorizontalProbingCells = + (unsigned)std::ceil( _nProbingCellsLineGreen / 2 * ( horizontalLineWidth / ( horizontalLineWidth + verticalLineWidth ) ) ); + unsigned nVerticalProbingCells = _nProbingCellsLineGreen - nHorizontalProbingCells; + + _probingCellsLineGreen = std::vector( _nProbingCellsLineGreen ); + + // printf( "here" ); + + // Sample points on each side of the rectangle + std::vector side1 = linspace2D( _cornerUpperLeftGreen, _cornerLowerLeftGreen, nVerticalProbingCells ); + std::vector side2 = linspace2D( _cornerLowerLeftGreen, _cornerLowerRightGreen, nHorizontalProbingCells ); + std::vector side3 = linspace2D( _cornerLowerRightGreen, _cornerUpperRightGreen, nVerticalProbingCells ); + std::vector side4 = linspace2D( _cornerUpperRightGreen, _cornerUpperLeftGreen, nHorizontalProbingCells ); + + // printf( "here" ); + // Combine the points from each side + _probingCellsLineGreen.insert( _probingCellsLineGreen.end(), side1.begin(), side1.end() ); + _probingCellsLineGreen.insert( _probingCellsLineGreen.end(), side2.begin(), side2.end() ); + _probingCellsLineGreen.insert( _probingCellsLineGreen.end(), side3.begin(), side3.end() ); + _probingCellsLineGreen.insert( _probingCellsLineGreen.end(), side4.begin(), side4.end() ); +} + +void SNSolverHPC::ComputeQOIsGreenProbingLine() { + std::cout << "ComputeQOIsGreenProbingLine" << std::endl; + double verticalLineWidth = std::abs( _cornerUpperLeftGreen[1] - _cornerLowerLeftGreen[1] ); + double horizontalLineWidth = std::abs( _cornerUpperLeftGreen[0] - _cornerUpperRightGreen[0] ); + + double dl = 2 * ( horizontalLineWidth + verticalLineWidth ) / ( (double)_nProbingCellsLineGreen ); + double area = dl * _thicknessGreen; + double a_g = 0; + double l_max = _nProbingCellsLineGreen * dl; + + for( unsigned i = 0; i < _nProbingCellsLineGreen; i++ ) { // Loop over probing cells + _absorptionValsIntegrated[i] = + ( _sigmaT[_probingCellsLineGreen[i]] - _sigmaS[_probingCellsLineGreen[i]] ) * _scalarFlux[_probingCellsLineGreen[i]] * area; + a_g += _absorptionValsIntegrated[i] / (double)_nProbingCellsLineGreen; + } + for( unsigned i = 0; i < _nProbingCellsLineGreen; i++ ) { // Loop over probing cells + _varAbsorptionValsIntegrated[i] = dl / l_max * ( a_g - _absorptionValsIntegrated[i] ) * ( a_g - _absorptionValsIntegrated[i] ); + } +} + +std::vector SNSolverHPC::linspace2D( const std::vector& start, const std::vector& end, unsigned num_points ) { + /** + * Generate a 2D linspace based on the start and end points with a specified number of points. + * + * @param start vector of starting x and y coordinates + * @param end vector of ending x and y coordinates + * @param num_points number of points to generate + * + * @return vector of unsigned integers representing the result + */ + + std::vector result; + result.resize( num_points ); + double stepX = ( end[0] - start[0] ) / ( num_points - 1 ); + double stepY = ( end[1] - start[1] ) / ( num_points - 1 ); + + for( unsigned i = 0; i < num_points; ++i ) { + double x = start[0] + i * stepX; + double y = start[1] + i * stepY; + + result[i] = _mesh->GetCellOfKoordinate( x, y ); + } + + return result; +} diff --git a/src/solvers/solverbase.cpp b/src/solvers/solverbase.cpp index 67f38e86..fbae12f6 100644 --- a/src/solvers/solverbase.cpp +++ b/src/solvers/solverbase.cpp @@ -37,7 +37,7 @@ SolverBase::SolverBase( Config* settings ) { // build slope related params _reconstructor = new Reconstructor( settings ); // Not used! - _reconsOrder = _reconstructor->GetReconsOrder(); + _reconsOrder = _reconstructor->GetSpatialOrder(); _interfaceMidPoints = _mesh->GetInterfaceMidPoints(); _cellMidPoints = _mesh->GetCellMidPoints(); @@ -61,6 +61,9 @@ SolverBase::SolverBase( Config* settings ) { if( _settings->GetIsCSD() ) { _nIter = _nEnergies - 1; // Since CSD does not go the last energy step } + else { + _nIter++; + } // setup problem and store frequently used params @@ -124,7 +127,7 @@ void SolverBase::Solve() { // Preprocessing before first pseudo time step SolverPreprocessing(); - unsigned rkStages = _settings->GetRKStages(); + unsigned rkStages = _settings->GetTemporalOrder(); // Create Backup solution for Runge Kutta VectorVector solRK0 = _sol; @@ -132,6 +135,10 @@ void SolverBase::Solve() { std::chrono::duration duration; // Loop over energies (pseudo-time of continuous slowing down approach) for( unsigned iter = 0; iter < _nIter; iter++ ) { + if( iter == _nIter - 1 ) { // last iteration + _dT = _settings->GetTEnd() - iter * _dT; + } + if( rkStages == 2 ) solRK0 = _sol; for( unsigned rkStep = 0; rkStep < rkStages; ++rkStep ) { @@ -175,7 +182,7 @@ void SolverBase::Solve() { DrawPostSolverOutput(); } -void SolverBase::RKUpdate( VectorVector sol_0, VectorVector sol_rk ) { +void SolverBase::RKUpdate( VectorVector& sol_0, VectorVector& sol_rk ) { #pragma omp parallel for for( unsigned i = 0; i < _nCells; ++i ) { _sol[i] = 0.5 * ( sol_0[i] + sol_rk[i] ); diff --git a/src/toolboxes/reconstructor.cpp b/src/toolboxes/reconstructor.cpp index 196e3443..f1d556b8 100644 --- a/src/toolboxes/reconstructor.cpp +++ b/src/toolboxes/reconstructor.cpp @@ -1,7 +1,7 @@ #include "toolboxes/reconstructor.hpp" #include "common/config.hpp" -Reconstructor::Reconstructor( Config* settings ) { _reconsOrder = settings->GetReconsOrder(); } +Reconstructor::Reconstructor( Config* settings ) { _reconsOrder = settings->GetSpatialOrder(); } double FortSign( double a, double b ) { if( b > 0.0 ) return std::fabs( a );