Organized project properly and added some property testing

.gitignore

@@ -77,5 +77,8 @@ Manifest.toml
 notes/*
 # My terrible, newbie attempt at version control (yes, shouldve jsut learned git)
 backups/*
-# Directory for holding the code for various bhmie implementations
+# Directory for caching test data
-.bhmielibs/*
+*/.cache/*
+# For now, performance will be tracked in a benchmarks folder
+*benchmarks
+*benchmarks/*

Project.toml

@@ -1,6 +1,9 @@
+name = "nanoconc"
+uuid = "9a947172-b1ea-4b16-84a6-f3d50752424d"
+authors = ["Cian Hughes <chughes000@gmail.com>"]
+version = "0.1.0"
+
 [deps]
-AirspeedVelocity = "1c8270ee-6884-45cc-9545-60fa71ec23e4"
-BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
 CSV = "336ed68f-0bac-5ca0-87d4-7b16caf5d00b"
 DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
 Debugger = "31a5f54b-26ea-5ae9-a837-f05ce5417438"

anchors.jl

@@ -0,0 +1,7 @@
+module Anchors
+
+const ROOT_DIR = @__DIR__
+const SRC_DIR = joinpath(ROOT_DIR, "src/")
+const TEST_DIR = joinpath(ROOT_DIR, "test/")
+
+end

test/Project.toml

@@ -0,0 +1,5 @@
+[deps]
+AirspeedVelocity = "1c8270ee-6884-45cc-9545-60fa71ec23e4"
+BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
+PropCheck = "ca382230-33be-11e9-0059-d981d03070e4"
+Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

test/benchmarks.jl

@@ -0,0 +1,83 @@
+module Benchmarks
+
+include("../anchors.jl")
+include("ffi_wraps.jl")
+
+import .Anchors.ROOT_DIR
+import .FFIWraps: bhmie_c, bhmie_fortran, bhmie_fortran77
+using BenchmarkTools
+
+include("$ROOT_DIR/src/miemfp.jl")
+
+function bench_vs_ffi()
+    # Fixed testing values
+    nang = UInt32(2)  # Example number of angles
+
+    c_result = @benchmark bhmie_c(x, cxref, nang, cxs1, cxs2) setup=(
+        x = rand(Float32);
+        cxref = rand(ComplexF32);
+        nang = UInt32($nang);
+        cxs1 = rand(ComplexF32, $nang);
+        cxs2 = rand(ComplexF32, $nang);
+    )
+
+    f_result = @benchmark bhmie_fortran(x, cxref, nang, cxs1, cxs2) setup=(
+        x = rand(Float32);
+        cxref = rand(ComplexF32);
+        nang = Int32($nang);
+        cxs1 = rand(ComplexF32, $nang);
+        cxs2 = rand(ComplexF32, $nang);
+    )
+
+    f77_result = @benchmark bhmie_fortran77(x, cxref, nang, cxs1, cxs2) setup=(
+        x = rand(Float32);
+        cxref = rand(ComplexF32);
+        nang = Int32($nang);
+        cxs1 = rand(ComplexF32, $nang);
+        cxs2 = rand(ComplexF32, $nang);
+    )
+
+    j_result = @benchmark miemfp.bhmie(Float64(x), ComplexF64(cxref), nang) setup=(
+        x = rand(Float32);
+        cxref = rand(ComplexF32);
+        nang = UInt32($nang);
+    )
+
+    println("\nC Implementation")
+    display(c_result)
+    println("\nFortran Implementation")
+    display(f_result)
+    println("\nFortran 77 Implementation")
+    display(f77_result)
+    println("\nJulia Implementation")
+    display(j_result)
+
+    return c_result, f_result, f77_result, j_result
+end
+
+end
+
+if abspath(PROGRAM_FILE) == @__FILE__
+    result = Benchmarks.bench_vs_ffi()
+
+    include("../anchors.jl")
+    import .Anchors.ROOT_DIR
+    using Pkg
+
+    current_package_version = Pkg.TOML.parsefile("$ROOT_DIR/Project.toml")["version"]
+
+    function display_to_file(io, x)
+        show(IOContext(io, :limit => false), "text/plain", x)
+    end
+
+    open("$ROOT_DIR/benchmarks/$current_package_version.txt", "w") do io
+        println(io, "C Implementation")
+        display_to_file(io, result[1])
+        println(io, "\n\nFortran Implementation")
+        display_to_file(io, result[2])
+        println(io, "\n\nFortran 77 Implementation")
+        display_to_file(io, result[3])
+        println(io, "\n\nJulia Implementation")
+        display_to_file(io, result[4])
+    end
+end

test/build_ffi.sh

@@ -10,9 +10,12 @@ base_url="http://scatterlib.wikidot.com/local--files/codes/"
 # Next, let's create a list of the codebases we want to pull
 codebases=("bhmie-f.zip" "bhmie-c.zip")
 
-# Then, let's create a directory to pull the codebases to
+# if bhmie_dir containers bhmie-c/bhmie.so, bhmie-f/bhmie.so, and bhmie-f/bhmie_f77.so then we can skip the build
-bhmie_dir=".bhmielibs"
+bhmie_dir=$1
-mkdir -p $bhmie_dir
+if [ -f $bhmie_dir/bhmie-c/bhmie.so ] && [ -f $bhmie_dir/bhmie-f/bhmie.so ] && [ -f $bhmie_dir/bhmie-f/bhmie_f77.so ]; then
+    echo "bhmie-c/bhmie.so, bhmie-f/bhmie.so, and bhmie-f/bhmie_f77.so already exist. Skipping build."
+    exit 0
+fi
 
 # Now, let's pull the codebases
 for codebase in ${codebases[@]}; do

test/ffi_wraps.jl

@@ -1,10 +1,21 @@
-run(`scripts/build_other_impls.sh`)
+module FFIWraps
 
-# using Debugger
+include("../anchors.jl")
-using BenchmarkTools
 
-include("src/miemfp.jl")
+if !@isdefined TEST_DIR
-# using miemfp
+    include("../anchors.jl")
+    import .Anchors: TEST_DIR
+end
+
+BHMIELIBS_DIR = joinpath(TEST_DIR, ".cache/bhmie-libs/")
+
+function __init__()
+    build_script = joinpath(TEST_DIR, "build_ffi.sh")
+
+    mkpath(BHMIELIBS_DIR)
+
+    run(`$build_script $BHMIELIBS_DIR`)
+end
 
 function bhmie_c(x::Float32, cxref::ComplexF32, nang::UInt32, cxs1::Vector{ComplexF32}, cxs2::Vector{ComplexF32})
     # Pre-allocate memory for the output variables
@@ -17,7 +28,7 @@ function bhmie_c(x::Float32, cxref::ComplexF32, nang::UInt32, cxs1::Vector{Compl
     # For example, if they need to be of size `nang`, you should verify or resize them accordingly
 
     # Call the C function
-    ccall((:bhmie, ".bhmielibs/bhmie-c/bhmie.so"), Cvoid,
+    ccall((:bhmie, "$BHMIELIBS_DIR/bhmie-c/bhmie.so"), Cvoid,
         (Float32, ComplexF32, UInt32, Vector{ComplexF32}, Vector{ComplexF32}, Ref{Float32}, Ref{Float32}, Ref{Float32}, Ref{Float32}),
         x, cxref, nang, cxs1, cxs2, qext, qsca, qback, gsca)
 
@@ -33,7 +44,7 @@ function bhmie_fortran(x::Float32, refrel::ComplexF32, nang::Int32, s1::Vector{C
     gsca = Ref{Float32}(0.0)
 
     # Call the Fortran subroutine
-    ccall((:bhmie_, ".bhmielibs/bhmie-f/bhmie.so"), Cvoid,
+    ccall((:bhmie_, "$BHMIELIBS_DIR/bhmie-f/bhmie.so"), Cvoid,
           (Ref{Float32}, Ref{ComplexF32}, Ref{Int32}, Ptr{ComplexF32}, Ptr{ComplexF32}, Ref{Float32}, Ref{Float32}, Ref{Float32}, Ref{Float32}),
           x, refrel, nang, s1, s2, qext, qsca, qback, gsca)
 
@@ -49,7 +60,7 @@ function bhmie_fortran77(x::Float32, refrel::ComplexF32, nang::Int32, s1::Vector
     gsca = Ref{Float32}(0.0)
 
     # Call the Fortran subroutine
-    ccall((:bhmie_, ".bhmielibs/bhmie-f/bhmie.so"), Cvoid,
+    ccall((:bhmie_, "$BHMIELIBS_DIR/bhmie-f/bhmie.so"), Cvoid,
           (Ref{Float32}, Ref{ComplexF32}, Ref{Int32}, Ptr{ComplexF32}, Ptr{ComplexF32}, Ref{Float32}, Ref{Float32}, Ref{Float32}, Ref{Float32}),
           x, refrel, nang, s1, s2, qext, qsca, qback, gsca)
 
@@ -57,44 +68,4 @@ function bhmie_fortran77(x::Float32, refrel::ComplexF32, nang::Int32, s1::Vector
     return qext[], qsca[], qback[], gsca[]
 end
 
-# Fixed testing values
+end
-nang = UInt32(2)  # Example number of angles
-
-c_result = @benchmark bhmie_c(x, cxref, nang, cxs1, cxs2) setup=(
-    x = rand(Float32);
-    cxref = rand(ComplexF32);
-    nang = UInt32($nang);
-    cxs1 = rand(ComplexF32, $nang);
-    cxs2 = rand(ComplexF32, $nang);
-)
-
-f_result = @benchmark bhmie_fortran(x, cxref, nang, cxs1, cxs2) setup=(
-    x = rand(Float32);
-    cxref = rand(ComplexF32);
-    nang = Int32($nang);
-    cxs1 = rand(ComplexF32, $nang);
-    cxs2 = rand(ComplexF32, $nang);
-)
-
-f77_result = @benchmark bhmie_fortran77(x, cxref, nang, cxs1, cxs2) setup=(
-    x = rand(Float32);
-    cxref = rand(ComplexF32);
-    nang = Int32($nang);
-    cxs1 = rand(ComplexF32, $nang);
-    cxs2 = rand(ComplexF32, $nang);
-)
-
-j_result = @benchmark miemfp.bhmie(Float64(x), ComplexF64(cxref), nang) setup=(
-    x = rand(Float32);
-    cxref = rand(ComplexF32);
-    nang = UInt32($nang);
-)
-
-println("\nC Implementation")
-display(c_result)
-println("\nFortran Implementation")
-display(f_result)
-println("\nFortran 77 Implementation")
-display(f77_result)
-println("\nJulia Implementation")
-display(j_result)

test/miemfp_tests.jl

@@ -1,31 +1,99 @@
 using Test
+using PropCheck
+
+include("../anchors.jl")
+
+import .Anchors: TEST_DIR, SRC_DIR
 
 if !@isdefined TestUtils
-    include("testutils.jl")
+    include(joinpath(TEST_DIR, "testutils.jl"))
 end
-if !@isdefined nanoconc
+if !@isdefined miemfp
-    include("../src/nanoconc.jl")
+    include(joinpath(SRC_DIR, "miemfp.jl"))
+end
+if !@isdefined FFIWraps
+    include(joinpath(TEST_DIR, "ffi_wraps.jl"))
 end
 
-@testset "miemfp" begin
+# function julia_vs_c(args::Tuple{Float64, Float64, Float64, UInt32, Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}})
+#     x, cxref_re, cxref_im, nang, cxs1_re, cxs1_im, cxs2_re, cxs2_im = args
+function julia_vs_c(x, cxref_re, cxref_im, nang, cxs1_re, cxs1_im, cxs2_re, cxs2_im)
+    cxref, cxs1, cxs2 = ComplexF64(cxref_re, cxref_im), ComplexF64.(cxs1_re, cxs1_im), ComplexF64.(cxs2_re, cxs2_im)
+    x_c, cxref_c, nang_c, cxs1_c, cxs2_c = Float32(x), ComplexF32(cxref), UInt32(nang), ComplexF32.(cxs1), ComplexF32.(cxs2)
+    return isapprox(
+        miemfp.bhmie(x, cxref, nang),
+        FFIWraps.bhmie_c(x_c, cxref_c, nang_c, cxs1_c, cxs2_c),
+        rtol=0.1,
+    )
+end
+
+# function julia_vs_fortran(args::Tuple{Float64, Float64, Float64, UInt32, Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}})
+#     x, cxref_re, cxref_im, nang, cxs1_re, cxs1_im, cxs2_re, cxs2_im = args
+function julia_vs_fortran(x, cxref_re, cxref_im, nang, cxs1_re, cxs1_im, cxs2_re, cxs2_im)
+    cxref, cxs1, cxs2 = ComplexF64(cxref_re, cxref_im), ComplexF64.(cxs1_re, cxs1_im), ComplexF64.(cxs2_re, cxs2_im)
+    x_f, cxref_f, nang_f, cxs1_f, cxs2_f = Float32(x), ComplexF32(cxref), Int32(nang), ComplexF32.(cxs1), ComplexF32.(cxs2)
+    b = miemfp.bhmie(x, cxref, nang)
+    f = FFIWraps.bhmie_fortran(x_f, cxref_f, nang_f, cxs1_f, cxs2_f)
+    # open("bhmie_julia_vs_fortran.txt", "a") do io
+    #     println(io, "julia: ", b)
+    #     println(io, "fortran: ", f)
+    # end
+    # return is_approx(b, f)
+    return isapprox(
+        miemfp.bhmie(x, cxref, nang),
+        FFIWraps.bhmie_fortran(x_f, cxref_f, nang_f, cxs1_f, cxs2_f),
+        rtol=0.1,
+    )
+end
+
+# function julia_vs_fortran77(args::Tuple{Float64, Float64, Float64, UInt32, Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}})
+#     x, cxref_re, cxref_im, nang, cxs1_re, cxs1_im, cxs2_re, cxs2_im = args
+function julia_vs_fortran77(x, cxref_re, cxref_im, nang, cxs1_re, cxs1_im, cxs2_re, cxs2_im)
+    cxref, cxs1, cxs2 = ComplexF64(cxref_re, cxref_im), ComplexF64.(cxs1_re, cxs1_im), ComplexF64.(cxs2_re, cxs2_im)
+    x_f, cxref_f, nang_f, cxs1_f, cxs2_f = Float32(x), ComplexF32(cxref), Int32(nang), ComplexF32.(cxs1), ComplexF32.(cxs2)
+    return isapprox(
+        miemfp.bhmie(x, cxref, nang),
+        FFIWraps.bhmie_fortran77(x_f, cxref_f, nang_f, cxs1_f, cxs2_f),
+        rtol=0.1,
+    )
+end
+
+f64_gen = PropCheck.itype(Float64)
+UInt32_gen = PropCheck.itype(UInt32)
+f64_vec_gen = PropCheck.vector(isample(1:100), f64_gen)
+bhmie_gen = PropCheck.interleave(
+    f64_gen,
+    f64_gen,
+    f64_gen,
+    UInt32_gen,
+    f64_vec_gen,
+    f64_vec_gen,
+    f64_vec_gen,
+    f64_vec_gen,
+)
+
+@testset "bhmie" begin
     @testset "miemfp.bhmie" begin
-        TestUtils.test_from_serialized(
+        c_check = PropCheck.check(julia_vs_c, bhmie_gen)
-            nanoconc.miemfp.bhmie,
+        c_result = c_check == true
-            "test/data/Main.nanoconc.miemfp.bhmie.ser"
+        if !c_result
-        )
+            println("Fail vs C, PropCheck:")
-    end
+            display(c_check)
-
+        end
-    @testset "miemfp.mfp" begin
+        @test c_result
-        TestUtils.test_from_serialized(
+        f_check = PropCheck.check(julia_vs_fortran, bhmie_gen)
-            nanoconc.miemfp.mfp,
+        f_result = f_check == true
-            "test/data/Main.nanoconc.miemfp.mfp.ser"
+        if !f_result
-        )
+            println("Fail vs Fortran, PropCheck:")
-    end
+            display(f_check)
-
+        end
-    @testset "miemfp.qbare" begin
+        @test f_result
-        TestUtils.test_from_serialized(
+        f77_check = PropCheck.check(julia_vs_fortran77, bhmie_gen)
-            nanoconc.miemfp.qbare,
+        f77_result = f77_check == true
-            "test/data/Main.nanoconc.miemfp.qbare.ser"
+        if !f77_result
-        )
+            println("Fail vs Fortran77, PropCheck:")
+            display(f77_check)
+        end
+        @test f77_result
     end
 end

test/runtests.jl

@@ -3,6 +3,9 @@ using Test
 include("testutils.jl")
 include("../src/nanoconc.jl")
 
-include("nanoconc_tests.jl")
+# include("nanoconc_tests.jl")
 include("miemfp_tests.jl")
-include("quantumcalc_tests.jl")
+# include("quantumcalc_tests.jl")
+include("benchmarks.jl")
+
+Benchmarks.bench_vs_ffi()