Notebook

Start-to-Finish Example: Unit Testing `GiRaFFE_NRPy`: Boundary Conditions¶

Author: Patrick Nelson¶

Edits by Terrence Pierre Jacques¶

This module Validates the Boundary Conditions routines for `GiRaFFE_NRPy`.¶

Notebook Status: Validated

Validation Notes: This module will validate the routines in Tutorial-GiRaFFE_NRPy-BCs.

NRPy+ Source Code for this module:¶

GiRaFFE_NRPy/GiRaFFE_NRPy_BCs.py [tutorial] Generates the driver apply boundary conditions to the vector potential and velocity.

Introduction:¶

This notebook validates the C code used to apply boundary conditions to components of the vector potential and valencia velocity.

It is, in general, good coding practice to unit test functions individually to verify that they produce the expected and intended output. We will generate test data with arbitrarily-chosen analytic functions and calculate gridfunctions at the cell centers on a small numeric grid. We will then compute the values for the ghost zones in two ways: first with the boundary condition C code driver, then we compute them analytically.

When this notebook is run, the significant digits of agreement between the approximate and exact values in the ghost zones will be evaluated. If the agreement falls below a thresold, the point, quantity, and level of agreement are reported here.

Table of Contents¶

$$\label{toc}$$

This notebook is organized as follows

Step 1: Set up core functions and parameters for unit testing the BCs algorithm
1. Step 1.a Write expressions for the gridfunctions we will test
2. Step 1.b Generate C functions to calculate the gridfunctions
3. Step 1.c Set free parameters in the code
Step 2: BCs_unit_test.c: The Main C Code
1. Step 2.a: Compile and run the code
Step 3: Code validation: Verify that relative error in numerical solution converges to zero at the expected order
Step 4: Output this notebook to $\LaTeX$-formatted PDF file

Step 1: Set up core functions and parameters for unit testing the BCs algorithm [Back to top]¶

$$\label{setup}$$

We'll start by appending the relevant paths to sys.path so that we can access sympy modules in other places. Then, we'll import NRPy+ core functionality and set up a directory in which to carry out our test.

In [1]:

import os, sys, shutil         # Standard Python modules for multiplatform OS-level functions
# First, we'll add the parent directory to the list of directories Python will check for modules.
nrpy_dir_path = os.path.join("..")
if nrpy_dir_path not in sys.path:
    sys.path.append(nrpy_dir_path)
nrpy_dir_path = os.path.join("..","..")
if nrpy_dir_path not in sys.path:
    sys.path.append(nrpy_dir_path)

from outputC import outCfunction, lhrh # NRPy+: Core C code output module
import sympy as sp               # SymPy: The Python computer algebra package upon which NRPy+ depends
import finite_difference as fin  # NRPy+: Finite difference C code generation module
import NRPy_param_funcs as par   # NRPy+: Parameter interface
import grid as gri               # NRPy+: Functions having to do with numerical grids
import indexedexp as ixp         # NRPy+: Symbolic indexed expression (e.g., tensors, vectors, etc.) support
import reference_metric as rfm   # NRPy+: Reference metric support
import cmdline_helper as cmd     # NRPy+: Multi-platform Python command-line interface

Ccodesdir = "Start-to-Finish-UnitTests/BCs_UnitTest/"

# First remove C code output directory if it exists
# Courtesy https://stackoverflow.com/questions/303200/how-do-i-remove-delete-a-folder-that-is-not-empty
shutil.rmtree(Ccodesdir, ignore_errors=True)
# Then create a fresh directory
cmd.mkdir(Ccodesdir)

outdir = os.path.join(Ccodesdir,"output/")
cmd.mkdir(outdir)

thismodule = "Start_to_Finish_UnitTest-GiRaFFE_NRPy-BCs"

# Set the finite-differencing order to 2
par.set_parval_from_str("finite_difference::FD_CENTDERIVS_ORDER", 2)

Step 1.a: Write expressions for the gridfunctions we will test [Back to top]¶

$$\label{expressions}$$

Now, we'll choose some functions with arbitrary forms to generate test data. We'll need to set seven gridfunctions, so expressions are being pulled from several previously written unit tests.

\begin{align} A_x &= dy + ez + f \\ A_y &= mx + nz + p \\ A_z &= sx + ty + u. \\ \bar{v}^x &= ax + by + cz \\ \bar{v}^y &= bx + cy + az \\ \bar{v}^z &= cx + ay + bz \\ [\sqrt{\gamma} \Phi] &= 1 - (x+2y+z) \\ \end{align}

In [2]:

a,b,c,d,e,f,g,h,l,m,n,o,p,q,r,s,t,u = par.Cparameters("REAL",thismodule,["a","b","c","d","e","f","g","h","l","m","n","o","p","q","r","s","t","u"],1e300)
M_PI  = par.Cparameters("#define",thismodule,["M_PI"], "")

AD = ixp.register_gridfunctions_for_single_rank1("EVOL","AD",DIM=3)
ValenciavU = ixp.register_gridfunctions_for_single_rank1("AUXEVOL","ValenciavU",DIM=3)
psi6Phi = gri.register_gridfunctions("EVOL","psi6Phi")

par.set_parval_from_str("reference_metric::CoordSystem","Cartesian")
rfm.reference_metric()
x = rfm.xx_to_Cart[0]
y = rfm.xx_to_Cart[1]
z = rfm.xx_to_Cart[2]

AD[0] = d*y + e*z + f
AD[1] = m*x + n*z + o
AD[2] = s*x + t*y + u

ValenciavU[0] = a#*x + b*y + c*z
ValenciavU[1] = b#*x + c*y + a*z
ValenciavU[2] = c#*x + a*y + b*z

psi6Phi = sp.sympify(1) - (x + sp.sympify(2)*y + z)

Step 1.b: Generate C functions to calculate the gridfunctions [Back to top]¶

$$\label{ccodekernels}$$

Here, we will use the NRPy+ function outCfunction() to generate C code that will calculate our metric gridfunctions over an entire grid; note that we call the function twice, once over just the interior points, and once over all points. This will allow us to compare against exact values in the ghostzones. We will also call the function to generate the boundary conditions function we are testing.

In [3]:

metric_gfs_to_print = [\
                       lhrh(lhs=gri.gfaccess("evol_gfs","AD0"),rhs=AD[0]),\
                       lhrh(lhs=gri.gfaccess("evol_gfs","AD1"),rhs=AD[1]),\
                       lhrh(lhs=gri.gfaccess("evol_gfs","AD2"),rhs=AD[2]),\
                       lhrh(lhs=gri.gfaccess("auxevol_gfs","ValenciavU0"),rhs=ValenciavU[0]),\
                       lhrh(lhs=gri.gfaccess("auxevol_gfs","ValenciavU1"),rhs=ValenciavU[1]),\
                       lhrh(lhs=gri.gfaccess("auxevol_gfs","ValenciavU2"),rhs=ValenciavU[2]),\
                       lhrh(lhs=gri.gfaccess("evol_gfs","psi6Phi"),rhs=psi6Phi),\
                      ]

desc = "Calculate test data on the interior grid for boundary conditions"
name = "calculate_test_data"
outCfunction(
    outfile  = os.path.join(Ccodesdir,name+".h"), desc=desc, name=name,
    params   ="const paramstruct *restrict params,REAL *restrict xx[3],REAL *restrict auxevol_gfs,REAL *restrict evol_gfs",
    body     = fin.FD_outputC("returnstring",metric_gfs_to_print,params="outCverbose=False"),
    loopopts="InteriorPoints,Read_xxs")

desc = "Calculate test data at all points for comparison"
name = "calculate_test_data_exact"
outCfunction(
    outfile  = os.path.join(Ccodesdir,name+".h"), desc=desc, name=name,
    params   ="const paramstruct *restrict params,REAL *restrict xx[3],REAL *restrict auxevol_gfs,REAL *restrict evol_gfs",
    body     = fin.FD_outputC("returnstring",metric_gfs_to_print,params="outCverbose=False"),
    loopopts="AllPoints,Read_xxs")

import GiRaFFE_NRPy.GiRaFFE_NRPy_BCs as BC
BC.GiRaFFE_NRPy_BCs(os.path.join(Ccodesdir,"boundary_conditions"))

Output C function calculate_test_data() to file Start-to-Finish-UnitTests/BCs_UnitTest/calculate_test_data.h
Output C function calculate_test_data_exact() to file Start-to-Finish-UnitTests/BCs_UnitTest/calculate_test_data_exact.h

Step 1.c: Set free parameters in the code [Back to top]¶

$$\label{free_parameters}$$

We also need to create the files that interact with NRPy's C parameter interface.

In [4]:

# Step 3.d.i: Generate declare_Cparameters_struct.h, set_Cparameters_default.h, and set_Cparameters[-SIMD].h
# par.generate_Cparameters_Ccodes(os.path.join(out_dir))

# Step 3.d.ii: Set free_parameters.h
with open(os.path.join(Ccodesdir,"free_parameters.h"),"w") as file:
    file.write("""
// Override parameter defaults with values based on command line arguments and NGHOSTS.
params.Nxx0 = atoi(argv[1]);
params.Nxx1 = atoi(argv[2]);
params.Nxx2 = atoi(argv[3]);
params.Nxx_plus_2NGHOSTS0 = params.Nxx0 + 2*NGHOSTS;
params.Nxx_plus_2NGHOSTS1 = params.Nxx1 + 2*NGHOSTS;
params.Nxx_plus_2NGHOSTS2 = params.Nxx2 + 2*NGHOSTS;
// Step 0d: Set up space and time coordinates
// Step 0d.i: Declare \Delta x^i=dxx{0,1,2} and invdxx{0,1,2}, as well as xxmin[3] and xxmax[3]:
const REAL xxmin[3] = {-1.0,-1.0,-1.0};
const REAL xxmax[3] = { 1.0, 1.0, 1.0};

params.dxx0 = (xxmax[0] - xxmin[0]) / ((REAL)params.Nxx_plus_2NGHOSTS0-1.0);
params.dxx1 = (xxmax[1] - xxmin[1]) / ((REAL)params.Nxx_plus_2NGHOSTS1-1.0);
params.dxx2 = (xxmax[2] - xxmin[2]) / ((REAL)params.Nxx_plus_2NGHOSTS2-1.0);
params.invdx0 = 1.0 / params.dxx0;
params.invdx1 = 1.0 / params.dxx1;
params.invdx2 = 1.0 / params.dxx2;
\n""")

# Generates declare_Cparameters_struct.h, set_Cparameters_default.h, and set_Cparameters[-SIMD].h
par.generate_Cparameters_Ccodes(os.path.join(Ccodesdir))

Step 2: `BCs_unit_test.c`: The Main C Code [Back to top]¶

$$\label{mainc}$$

Here we compare the results of our boundary conditions, apply_bcs_potential() and apply_bcs_velocity(), against the exact results, looping over the entire numerical grid.

In [5]:

%%writefile $Ccodesdir/BCs_unit_test.c

// These are common packages that we are likely to need.
#include "stdio.h"
#include "stdlib.h"
#include "math.h"
#include "string.h" // Needed for strncmp, etc.
#include "stdint.h" // Needed for Windows GCC 6.x compatibility
#include <time.h>   // Needed to set a random seed.

#define REAL double
#include "declare_Cparameters_struct.h"

const int NGHOSTS = 3;

REAL a,b,c,d,e,f,g,h,l,m,n,o,p,q,r,s,t,u;

// Standard NRPy+ memory access:
#define IDX4S(g,i,j,k) \
( (i) + Nxx_plus_2NGHOSTS0 * ( (j) + Nxx_plus_2NGHOSTS1 * ( (k) + Nxx_plus_2NGHOSTS2 * (g) ) ) )

// Standard formula to calculate significant digits of agreement:
#define SDA(a,b) 1.0-log10(2.0*fabs(a-b)/(fabs(a)+fabs(b)))

// Give gridfunctions their names:
#define VALENCIAVU0GF 0
#define VALENCIAVU1GF 1
#define VALENCIAVU2GF 2
#define NUM_AUXEVOL_GFS 3

#define AD0GF 0
#define AD1GF 1
#define AD2GF 2
#define STILDED0GF 3
#define STILDED1GF 4
#define STILDED2GF 5
#define PSI6PHIGF 6
#define NUM_EVOL_GFS 7

#include "calculate_test_data.h"
#include "calculate_test_data_exact.h"
#include "boundary_conditions/GiRaFFE_boundary_conditions.h"

int main(int argc, const char *argv[]) {
    paramstruct params;
#include "set_Cparameters_default.h"

    // Step 0c: Set free parameters, overwriting Cparameters defaults
    //          by hand or with command-line input, as desired.
#include "free_parameters.h"
#include "set_Cparameters-nopointer.h"

    // We'll define our grid slightly different from how we normally would. We let our outermost
    // ghostzones coincide with xxmin and xxmax instead of the interior of the grid. This means
    // that the ghostzone points will have identical positions so we can do convergence tests of them.
    // Step 0e: Set up cell-centered Cartesian coordinate grids
    REAL *xx[3];
    xx[0] = (REAL *)malloc(sizeof(REAL)*Nxx_plus_2NGHOSTS0);
    xx[1] = (REAL *)malloc(sizeof(REAL)*Nxx_plus_2NGHOSTS1);
    xx[2] = (REAL *)malloc(sizeof(REAL)*Nxx_plus_2NGHOSTS2);
    for(int j=0;j<Nxx_plus_2NGHOSTS0;j++) xx[0][j] = xxmin[0] + ((REAL)(j))*dxx0;
    for(int j=0;j<Nxx_plus_2NGHOSTS1;j++) xx[1][j] = xxmin[1] + ((REAL)(j))*dxx1;
    for(int j=0;j<Nxx_plus_2NGHOSTS2;j++) xx[2][j] = xxmin[2] + ((REAL)(j))*dxx2;

    //for(int i=0;i<Nxx_plus_2NGHOSTS0;i++) printf("xx[0][%d] = %.15e\n",i,xx[0][i]);

    // This is the array to which we'll write the NRPy+ variables.
    REAL *auxevol_gfs  = (REAL *)malloc(sizeof(REAL) * NUM_AUXEVOL_GFS * Nxx_plus_2NGHOSTS2 * Nxx_plus_2NGHOSTS1 * Nxx_plus_2NGHOSTS0);
    REAL *evol_gfs  = (REAL *)malloc(sizeof(REAL) * NUM_EVOL_GFS * Nxx_plus_2NGHOSTS2 * Nxx_plus_2NGHOSTS1 * Nxx_plus_2NGHOSTS0);
    // And another for exact data:
    REAL *auxevol_exact_gfs  = (REAL *)malloc(sizeof(REAL) * NUM_AUXEVOL_GFS * Nxx_plus_2NGHOSTS2 * Nxx_plus_2NGHOSTS1 * Nxx_plus_2NGHOSTS0);
    REAL *evol_exact_gfs  = (REAL *)malloc(sizeof(REAL) * NUM_EVOL_GFS * Nxx_plus_2NGHOSTS2 * Nxx_plus_2NGHOSTS1 * Nxx_plus_2NGHOSTS0);

    // Generate some random coefficients. Leave the random seed on its default for consistency between trials.
    a = (double)(rand()%20)/5.0;
    f = (double)(rand()%20)/5.0;
    m = (double)(rand()%20)/5.0;
    b = (double)(rand()%10-5)/100.0;
    c = (double)(rand()%10-5)/100.0;
    d = (double)(rand()%10-5)/100.0;
    g = (double)(rand()%10-5)/100.0;
    h = (double)(rand()%10-5)/100.0;
    l = (double)(rand()%10-5)/100.0;
    n = (double)(rand()%10-5)/100.0;
    o = (double)(rand()%10-5)/100.0;
    p = (double)(rand()%10-5)/100.0;

    // First, calculate the test data on our grid, along with the comparison:
    calculate_test_data(&params,xx,auxevol_gfs,evol_gfs);
    calculate_test_data_exact(&params,xx,auxevol_exact_gfs,evol_exact_gfs);

    // Run the BCs driver on the test data to fill in the ghost zones:
    apply_bcs_potential(&params,evol_gfs);
    apply_bcs_velocity(&params,auxevol_gfs);

    /*char filename[100];
    sprintf(filename,"out%d-numer.txt",Nxx0);
    FILE *out2D = fopen(filename, "w");
    for(int i2=0;i2<Nxx_plus_2NGHOSTS2;i2++) for(int i1=0;i1<Nxx_plus_2NGHOSTS1;i1++) for(int i0=0;i0<Nxx_plus_2NGHOSTS0;i0++) {
        // We print the difference between approximate and exact numbers.
        fprintf(out2D,"%.16e\t %e %e %e\n",
                //auxevol_gfs[IDX4S(VALENCIAVU2GF,i0,i1,i2)]-auxevol_exact_gfs[IDX4S(VALENCIAVU2GF,i0,i1,i2)],
                evol_gfs[IDX4S(AD2GF,i0,i1,i2)]-evol_exact_gfs[IDX4S(AD2GF,i0,i1,i2)],
                xx[0][i0],xx[1][i1],xx[2][i2]
                );
    }
    fclose(out2D);*/

    int all_agree = 1;

    for(int i0=0;i0<Nxx_plus_2NGHOSTS0;i0++){
        for(int i1=0;i1<Nxx_plus_2NGHOSTS1;i1++){
            for(int i2=0;i2<Nxx_plus_2NGHOSTS2;i2++){
                if(SDA(evol_gfs[IDX4S(AD0GF, i0,i1,i2)],evol_exact_gfs[IDX4S(AD0GF, i0,i1,i2)])<10.0){
                    printf("Quantity AD0 only agrees with the original GiRaFFE to %.2f digits at i0,i1,i2=%d,%d,%d!\n",
                           SDA(evol_gfs[IDX4S(AD0GF, i0,i1,i2)],evol_exact_gfs[IDX4S(AD0GF, i0,i1,i2)]),i0,i1,i2);
                    all_agree=0;
                }
                if(SDA(evol_gfs[IDX4S(AD1GF, i0,i1,i2)],evol_exact_gfs[IDX4S(AD1GF, i0,i1,i2)])<10.0){
                    printf("Quantity AD1 only agrees with the original GiRaFFE to %.2f digits at i0,i1,i2=%d,%d,%d!\n",
                           SDA(evol_gfs[IDX4S(AD1GF, i0,i1,i2)],evol_exact_gfs[IDX4S(AD1GF, i0,i1,i2)]),i0,i1,i2);
                    all_agree=0;
                }
                if(SDA(evol_gfs[IDX4S(AD2GF, i0,i1,i2)],evol_exact_gfs[IDX4S(AD2GF, i0,i1,i2)])<10.0){
                    printf("Quantity AD2 only agrees with the original GiRaFFE to %.2f digits at i0,i1,i2=%d,%d,%d!\n",
                           SDA(evol_gfs[IDX4S(AD2GF, i0,i1,i2)],evol_exact_gfs[IDX4S(AD2GF, i0,i1,i2)]),i0,i1,i2);
                    all_agree=0;
                }
                if(SDA(auxevol_gfs[IDX4S(VALENCIAVU0GF, i0,i1,i2)],auxevol_exact_gfs[IDX4S(VALENCIAVU0GF, i0,i1,i2)])<10.0){
                    printf("Quantity ValenciavU0 only agrees with the original GiRaFFE to %.2f digits at i0,i1,i2=%d,%d,%d!\n",
                           SDA(auxevol_gfs[IDX4S(VALENCIAVU0GF, i0,i1,i2)],auxevol_exact_gfs[IDX4S(VALENCIAVU0GF, i0,i1,i2)]),i0,i1,i2);
                    all_agree=0;
                }
                if(SDA(auxevol_gfs[IDX4S(VALENCIAVU1GF, i0,i1,i2)],auxevol_exact_gfs[IDX4S(VALENCIAVU1GF, i0,i1,i2)])<10.0){
                    printf("Quantity ValenciavU1 only agrees with the original GiRaFFE to %.2f digits at i0,i1,i2=%d,%d,%d!\n",
                           SDA(auxevol_gfs[IDX4S(VALENCIAVU1GF, i0,i1,i2)],auxevol_exact_gfs[IDX4S(VALENCIAVU1GF, i0,i1,i2)]),i0,i1,i2);
                    all_agree=0;
                }
                if(SDA(auxevol_gfs[IDX4S(VALENCIAVU2GF, i0,i1,i2)],auxevol_exact_gfs[IDX4S(VALENCIAVU2GF, i0,i1,i2)])<10.0){
                    printf("Quantity ValenciavU2 only agrees with the original GiRaFFE to %.2f digits at i0,i1,i2=%d,%d,%d!\n",
                           SDA(auxevol_gfs[IDX4S(VALENCIAVU2GF, i0,i1,i2)],auxevol_exact_gfs[IDX4S(VALENCIAVU2GF, i0,i1,i2)]),i0,i1,i2);
                    all_agree=0;
                }
                if(SDA(evol_gfs[IDX4S(PSI6PHIGF, i0,i1,i2)],evol_exact_gfs[IDX4S(PSI6PHIGF, i0,i1,i2)])<10.0){
                    printf("psi6Phi = %.15e,%.15e\n",evol_gfs[IDX4S(PSI6PHIGF, i0,i1,i2)],evol_exact_gfs[IDX4S(PSI6PHIGF, i0,i1,i2)]);
                    //printf("Quantity psi6Phi only agrees with the original GiRaFFE to %.2f digits at i0,i1,i2=%d,%d,%d!\n",
                    //       SDA(evol_gfs[IDX4S(PSI6PHIGF, i0,i1,i2)],evol_exact_gfs[IDX4S(PSI6PHIGF, i0,i1,i2)]),i0,i1,i2);
                    all_agree=0;
                }
            }
        }
    }
    if(all_agree) printf("All quantities agree at all points!\n");

}

Writing Start-to-Finish-UnitTests/BCs_UnitTest//BCs_unit_test.c

Step 2.a: Compile and run the code [Back to top]¶

$$\label{compile_run}$$

Now that we have our file, we can compile it and run the executable.

In [6]:

import time
results_file = "out.txt"


print("Now compiling, should take ~2 seconds...\n")
start = time.time()
cmd.C_compile(os.path.join(Ccodesdir,"BCs_unit_test.c"), os.path.join(outdir,"BCs_unit_test"))
end = time.time()
print("Finished in "+str(end-start)+" seconds.\n\n")

# Change to output directory
os.chdir(outdir)

print("Now running...\n")
start = time.time()
cmd.Execute(os.path.join("BCs_unit_test"),"2 2 2",file_to_redirect_stdout=results_file)
# To do a convergence test, we'll also need a second grid with twice the resolution.
# cmd.Execute(os.path.join("Validation","BCs_unit_test"),"9 9 9",file_to_redirect_stdout=os.path.join(out_dir,results_file))
end = time.time()
print("Finished in "+str(end-start)+" seconds.\n\n")

Now compiling, should take ~2 seconds...

Compiling executable...
(EXEC): Executing `gcc -Ofast -fopenmp -march=native -funroll-loops Start-to-Finish-UnitTests/BCs_UnitTest/BCs_unit_test.c -o Start-to-Finish-UnitTests/BCs_UnitTest/output/BCs_unit_test -lm`...
(BENCH): Finished executing in 0.6085014343261719 seconds.
Finished compilation.
Finished in 0.6146230697631836 seconds.


Now running...

(EXEC): Executing `taskset -c 0,1,2,3,4,5 ./BCs_unit_test 2 2 2`...
(BENCH): Finished executing in 0.21192646026611328 seconds.
Finished in 0.2279961109161377 seconds.

Now, we will interpret our output and verify that we produced the correct results.

In [7]:

with open(results_file,"r") as file:
    output = file.readline()
    print(output)
    if output!="All quantities agree at all points!\n": # If this isn't the first line of this file, something went wrong!
        sys.exit(1)

All quantities agree at all points!

Step 4: Output this notebook to $\LaTeX$-formatted PDF file [Back to top]¶

$$\label{latex_pdf_output}$$

The following code cell converts this Jupyter notebook into a proper, clickable $\LaTeX$-formatted PDF file. After the cell is successfully run, the generated PDF may be found in the root NRPy+ tutorial directory, with filename Tutorial-Start_to_Finish-GiRaFFE_NRPy-Metric_Face_Values.pdf (Note that clicking on this link may not work; you may need to open the PDF file through another means.)

In [8]:

import cmdline_helper as cmd    # NRPy+: Multi-platform Python command-line interface

# Change to NRPy directory
os.chdir("../../../")

cmd.output_Jupyter_notebook_to_LaTeXed_PDF("Tutorial-Start_to_Finish_UnitTest-GiRaFFE_NRPy-BCs",location_of_template_file=os.path.join(".."))

Created Tutorial-Start_to_Finish_UnitTest-GiRaFFE_NRPy-BCs.tex, and
    compiled LaTeX file to PDF file Tutorial-Start_to_Finish_UnitTest-
    GiRaFFE_NRPy-BCs.pdf

Start-to-Finish Example: Unit Testing GiRaFFE_NRPy: Boundary Conditions¶