Cleaned up and added comparisons against established implementations

This commit is contained in:
Cian Hughes
2024-02-02 17:22:28 +00:00
parent cf286c675b
commit 6948890a28
64 changed files with 141 additions and 11564 deletions

2
.gitignore vendored
View File

@@ -77,3 +77,5 @@ 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
.bhmielibs/*

3
.vscode/settings.json vendored Normal file
View File

@@ -0,0 +1,3 @@
{
"GitHooks.hooksDirectory": "/home/cianh/Programming/Git_Projects/nanoconc/.git/hooks"
}

View File

@@ -1,24 +0,0 @@
include("nanoconc.jl")
using PyPlot
particledata = ([2.5 1.],
[5. 1.],
[10. 1.],
[20. 1.],
[30. 1.],
[40. 1.],
[50. 1.],
[60. 1.],
[70. 1.],
[80. 1.],
[90. 1.],
[100. 1.])
mat = nanoconc.loadmaterial("Gold",disp=false)
predict(x) = nanoconc.abspredict(1.333,250.,900.,0x00001000,x,mat,1000.,1.)
spectra = predict.(particledata)
for spectrum in spectra
plot(spectrum[:,1],spectrum[:,2])
end
show()

View File

@@ -1,39 +0,0 @@
println("\nImporting modules...")
include("nanoconc.jl")
using PyPlot
using JLD2
println("Loading parameters...")
particledata = ([2.5 1.],
[5. 1.],
[10. 1.],
[20. 1.],
[30. 1.],
[40. 1.],
[50. 1.],
[60. 1.],
[70. 1.],
[80. 1.],
[90. 1.],
[100. 1.])
mat = nanoconc.loadmaterial("Gold",disp=false)
predict(x) = nanoconc.qpredict(1.333,250.,900.,0x00001000,x,mat)
println("Calculating spectra...")
print("Calculation time: ")
@time spectra = predict.(particledata)
println("Saving output...")
println("Plotting spectra...")
for spectrum in spectra
plot(spectrum[:,1],spectrum[:,2])
end
println("Displaying spectra...")
println("Script complete!")
show()
println()

28
20nm.jl
View File

@@ -1,28 +0,0 @@
println("\nImporting modules...")
include("nanoconc.jl")
using PyPlot
println("Loading parameters...")
particledata = [20. 1.]
mat = nanoconc.loadmaterial("Gold",disp=false)
predict(x) = nanoconc.abspredict(1.333,450.,650.,0x00001000,x,mat,1000.,1.)
println("Calculating spectra...")
print("Calculation time: ")
@time spectrum = predict(particledata)
println("Finding lambda max...")
#println("\nDebug data:\nfindmax = $(findmax(spectrum[:,2]))\nmaxind = $(findmax(spectrum[:,2])[2])")
lmax = spectrum[findmax(spectrum[:,2])[2],1]
println("lambda max = $(lmax)")
println("Plotting spectra...")
plot(spectrum[:,1],spectrum[:,2])
println("Displaying spectra...")
show()
println("Script complete!")

4096
abs.csv

File diff suppressed because it is too large Load Diff

View File

@@ -1,37 +0,0 @@
println("\nImporting modules...")
include("src/nanoconc.jl")
using CSV
using DataFrames
using Profile, ProfileVega
println("Loading parameters...")
refmed = 1.333 # refractive index of the medium
wavel1, wavel2 = 250., 900. # wavelengths to calculate between
numval = 0x00001000 # number of values to calculate in spectrum
particledata = [2.5 1.;
5. 1.;
10. 1.;
20. 1.;
30. 1.;
40. 1.;
50. 1.;
60. 1.;
70. 1.;
80. 1.;
90. 1.;
100. 1.]
materialdata = nanoconc.loadmaterial("Gold",disp=false)
ppml = 10000. # particles per ml nanoconc.quantumcalc.numtomol(ppml) -> conc
d0 = 1. # path length
params = (refmed,wavel1,wavel2,numval,particledata,materialdata,ppml,d0)
spectrum = nanoconc.abspredict(params...)
println("Calculating absorbance spectrum...")
print("Calculation time: ")
@time spectrum = nanoconc.abspredict(params...)
CSV.write("abs.csv", DataFrame(spectrum, :auto), writeheader=false)
# @profview nanoconc.abspredict(params...)

80
bhmielibs_test.jl Normal file
View File

@@ -0,0 +1,80 @@
run(`scripts/build_other_impls.sh`)
include("src/miemfp.jl")
# using miemfp
function bhmie_c(x::Float32, cxref::ComplexF32, nang::UInt32, cxs1::Vector{ComplexF32}, cxs2::Vector{ComplexF32})
# Pre-allocate memory for the output variables
qext = Ref{Float32}(0.0)
qsca = Ref{Float32}(0.0)
qback = Ref{Float32}(0.0)
gsca = Ref{Float32}(0.0)
# Ensure cxs1 and cxs2 have proper sizes, as expected by the C function
# 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,
(Float32, ComplexF32, UInt32, Ptr{ComplexF32}, Ptr{ComplexF32}, Ref{Float32}, Ref{Float32}, Ref{Float32}, Ref{Float32}),
x, cxref, nang, cxs1, cxs2, qext, qsca, qback, gsca)
# Return the output variables
return qext[], qsca[], qback[], gsca[]
end
function bhmie_fortran(x::Float32, refrel::ComplexF32, nang::Int32, s1::Vector{ComplexF32}, s2::Vector{ComplexF32})
# Pre-allocate output variables
qext = Ref{Float32}(0.0)
qsca = Ref{Float32}(0.0)
qback = Ref{Float32}(0.0)
gsca = Ref{Float32}(0.0)
# Call the Fortran subroutine
ccall((:bhmie_, ".bhmielibs/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)
# Return the modified values
return qext[], qsca[], qback[], gsca[]
end
function bhmie_fortran77(x::Float32, refrel::ComplexF32, nang::Int32, s1::Vector{ComplexF32}, s2::Vector{ComplexF32})
# Pre-allocate output variables
qext = Ref{Float32}(0.0)
qsca = Ref{Float32}(0.0)
qback = Ref{Float32}(0.0)
gsca = Ref{Float32}(0.0)
# Call the Fortran subroutine
ccall((:bhmie_, ".bhmielibs/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)
# Return the modified values
return qext[], qsca[], qback[], gsca[]
end
# Create test data
x = Float32(1.0) # Example value for x
cxref = ComplexF32(1.5, 0.5) # Example complex refractive index
nang = UInt32(2) # Example number of angles
# Example arrays for scattering amplitudes, initialized with dummy complex values
cxs1 = [ComplexF32(0.1, 0.2) for _ in 1:nang]
cxs2 = [ComplexF32(0.3, 0.4) for _ in 1:nang]
# Test C wrapper
qext, qsca, qback, gsca = bhmie_c(x, cxref, nang, cxs1, cxs2)
println("bhmie_c output: qext = $qext, qsca = $qsca, qback = $qback, gsca = $gsca")
# Test Fortran wrapper
qext, qsca, qback, gsca = bhmie_fortran(x, cxref, Int32(nang), cxs1, cxs2)
println("bhmie_fortran output: qext = $qext, qsca = $qsca, qback = $qback, gsca = $gsca")
# Test Fortran77 wrapper
qext, qsca, qback, gsca = bhmie_fortran77(x, cxref, Int32(nang), cxs1, cxs2)
println("bhmie_fortran77 output: qext = $qext, qsca = $qsca, qback = $qback, gsca = $gsca")
# Test Julia wrapper
qext, qsca, qback, gsca = miemfp.bhmie(Float64(x), ComplexF64(cxref), nang)
println("bhmie_julia output: qext = $qext, qsca = $qsca, qback = $qback, gsca = $gsca")

3594
prof.ipynb

File diff suppressed because one or more lines are too long

File diff suppressed because it is too large Load Diff

55
scripts/build_other_impls.sh Executable file
View File

@@ -0,0 +1,55 @@
#!/bin/bash
# This bash script will pull and compile the other implementations of the
# bhmie algorithm (found at http://scatterlib.wikidot.com/mie)
# These implementations will then be tested and benchmarked against the local
# implementation to ensure that they are correct and (ideally) faster.
# First, let's create a base URL for the scatterlib wiki codebases
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
bhmie_dir=".bhmielibs"
mkdir -p $bhmie_dir
# Now, let's pull the codebases
for codebase in ${codebases[@]}; do
wget -O $bhmie_dir/$codebase $base_url/$codebase
if [[ $codebase == *.zip ]]; then
unzip $bhmie_dir/$codebase -d $bhmie_dir/$codebase
fi
done
# Then, let's extract any .zip files to a directory of the same name (without the .zip, obviously)
for codebase in ${codebases[@]}; do
if [[ $codebase == *.zip ]]; then
unzip $bhmie_dir/$codebase -d $bhmie_dir/${codebase%.zip}
fi
done
# If the folder already existed in the archive, move that folder one level up
# and remove the now empty folder
for codebase in ${codebases[@]}; do
if [[ $codebase == *.zip ]]; then
folder=${codebase%.zip}
if [ -d $bhmie_dir/$folder/$folder ]; then
mv $bhmie_dir/$folder/$folder/* $bhmie_dir/$folder
rmdir $bhmie_dir/$folder/$folder
fi
fi
done
# Next, if this has all succeeded we can delete the .zip files
for codebase in ${codebases[@]}; do
if [[ $codebase == *.zip ]]; then
rm $bhmie_dir/$codebase
fi
done
# And, finally, we can compile the C, and Fortran implementations
cd $bhmie_dir
gcc -shared -fPIC -o bhmie-c/bhmie.so bhmie-c/bhmie.c bhmie-c/complex.c bhmie-c/nrutil.c -lm -Wno-builtin-declaration-mismatch -Wno-implicit-function-declaration
gfortran -shared -fPIC -o bhmie-f/bhmie.so bhmie-f/bhmie.f
gfortran -shared -fPIC -o bhmie-f/bhmie_f77.so bhmie-f/bhmie_f77.f

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

View File

@@ -1 +0,0 @@
Cian Hughes,cianh,cian-worklaptop,04.07.2023 11:49,file:///home/cianh/.config/libreoffice/4;

View File

@@ -1 +0,0 @@
Cian Hughes,cianh,cian-worklaptop,04.07.2023 14:52,file:///home/cianh/.config/libreoffice/4;

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

File diff suppressed because one or more lines are too long