-
Notifications
You must be signed in to change notification settings - Fork 1
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
- Loading branch information
1 parent
4eac287
commit 3f61c88
Showing
10 changed files
with
567 additions
and
157 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,16 @@ | ||
## Quick use test | ||
|
||
An introductory non-trivial example is given in `examples/run_plasma.jl`. This script will recreate Figure E1 from [Fouvry & Prunet (2021).](https://ui.adsabs.harvard.edu/abs/2022MNRAS.509.2443F/abstract) | ||
|
||
If you installed the library using the first (global) install option, just download this example [file](https://github.com/JuliaStellarDynamics/FiniteHilbertTransform.jl/blob/main/examples/run_plasma.jl) from the github repository. | ||
|
||
Run the code with the following command: | ||
``` | ||
$ julia /path/to/run_plasma.jl | ||
``` | ||
|
||
This example will first install some required libraries (`Plots`, `ArgParse`). These installations might take a few minutes when first called. | ||
|
||
The resulting plot will be created in the same folder as the test code under the name `plasmademo.png`. | ||
|
||
 |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Large diffs are not rendered by default.
Oops, something went wrong.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -1,57 +1,49 @@ | ||
""" | ||
tabLeg!_DAMPED(omg::Complex{Float64}, struct_tabLeg::LegendreFHT) | ||
Fill in all the Legendre arrays for a DAMPED mode, i.e., Im[w] < 0.0. | ||
# Arguments | ||
- `omg::Complex{Float64}`: Complex frequency. | ||
- `struct_tabLeg::LegendreFHT`: LegendreFHT structure. | ||
""" | ||
tabLeg!_DAMPED | ||
Function to fill in all the Legendre arrays | ||
for a DAMPED mode, i.e. Im[w] < 0.0 | ||
""" | ||
function tabLeg!_DAMPED(omg::ComplexF64, | ||
struct_tabLeg::LegendreFHT) | ||
# Description | ||
`tabLeg!_DAMPED` fills in all the Legendre arrays for a DAMPED mode, where Im[w] < 0.0. It computes the Legendre coefficients D_k(w) and, if needed, Legendre polynomials P_k(w) for the given complex frequency. | ||
# the container for D_k(w) | ||
""" | ||
function tabLeg!_DAMPED(omg::Complex{Float64}, struct_tabLeg::LegendreFHT) | ||
# The container for D_k(w) | ||
tabDLeg = struct_tabLeg.tabDLeg | ||
|
||
# the container for Q_k(w) | ||
# The container for Q_k(w) | ||
tabQLeg = struct_tabLeg.tabQLeg | ||
|
||
# now, compute Q_k(w) | ||
|
||
# fill the initial values of Q_k(w) | ||
val_0_Q = log(1.0-omg) - log(-1.0-omg) | ||
val_1_Q = 2.0 + omg*val_0_Q | ||
# Compute Q_k(w) | ||
# Fill the initial values of Q_k(w) | ||
val_0_Q = log(1.0 - omg) - log(-1.0 - omg) | ||
val_1_Q = 2.0 + omg * val_0_Q | ||
# Compute the rest of the k values | ||
tabQLeg!(omg, val_0_Q, val_1_Q, tabQLeg) | ||
|
||
# compute the rest of the k values | ||
tabQLeg!(omg,val_0_Q,val_1_Q,tabQLeg) | ||
|
||
# compute the value of H[-1,1,Re[w]]. | ||
# Compute the value of H[-1,1,Re[w]] | ||
val_H = Heaviside(real(omg)) | ||
|
||
# if we are outside of the interval [-1,1], we do not need to compute P_k(w) (it =0) | ||
if (val_H == 0.0) | ||
|
||
# If we are outside of the interval [-1,1], we do not need to compute P_k(w) (it =0) | ||
if val_H == 0.0 | ||
# Computing D_k(w) | ||
for k=0:(struct_tabLeg.Ku-1) # Loop over the Legendre indices | ||
tabDLeg[k+1] = tabQLeg[k+1] # Filling in tabDLeg. ATTENTION, to the shift of the arrays. | ||
for k in 0:(struct_tabLeg.Ku-1) # Loop over the Legendre indices | ||
tabDLeg[k+1] = tabQLeg[k+1] # Filling in tabDLeg | ||
end | ||
|
||
# if we are within the interval [-1,1], we need to compute P_k(w) | ||
else | ||
|
||
# the container for P_k(w) | ||
# The container for P_k(w) | ||
tabPLeg = struct_tabLeg.tabPLeg | ||
|
||
# Initial values of P_k(w) | ||
val_0_P = 1.0 + 0.0*im # Initial value for k = 0. ATTENTION, we make sure that the quantity is seen as complex. | ||
val_1_P = omg # Initial value for k = 1. | ||
##### | ||
tabPLeg!(omg,val_0_P,val_1_P,struct_tabLeg.Ku,tabPLeg) # Computing the tabPLeg | ||
########## | ||
val_0_P = 1.0 + 0.0im # Initial value for k = 0 | ||
val_1_P = omg # Initial value for k = 1 | ||
tabPLeg!(omg, val_0_P, val_1_P, struct_tabLeg.Ku, tabPLeg) # Compute tabPLeg | ||
|
||
# Computing the D_k(w) | ||
##### | ||
for k=0:(struct_tabLeg.Ku-1) # Loop over the Legendre indices | ||
tabDLeg[k+1] = tabQLeg[k+1] + 2.0*im*pi*val_H*tabPLeg[k+1] # Filling in tabDLeg. ATTENTION, to the shift of the arrays | ||
for k in 0:(struct_tabLeg.Ku-1) # Loop over the Legendre indices | ||
tabDLeg[k+1] = tabQLeg[k+1] + 2.0im * π * val_H * tabPLeg[k+1] # Filling in tabDLeg | ||
end | ||
##### | ||
end | ||
end |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -1,75 +1,59 @@ | ||
""" | ||
tabLeg!_NEUTRAL(omg::Complex{Float64}, struct_tabLeg::LegendreFHT) | ||
Fill in all the Legendre arrays for a NEUTRAL mode, i.e., Im[w] = 0.0. | ||
""" | ||
Function to fill in all the Legendre arrays for an NEUTRAL mode, i.e. Im[w] = 0.0 | ||
""" | ||
function tabLeg!_NEUTRAL(omg::ComplexF64, | ||
struct_tabLeg::LegendreFHT) | ||
# Arguments | ||
- `omg::Complex{Float64}`: Complex frequency. | ||
- `struct_tabLeg::LegendreFHT`: LegendreFHT structure. | ||
# The integral is computed in the Cauchy sense, hence the absolute values. | ||
# Description | ||
`tabLeg!_NEUTRAL` fills in all the Legendre arrays for a NEUTRAL mode, where Im[w] = 0.0. It computes the Legendre coefficients D_k(w) and, if needed, Legendre polynomials P_k(w) for the given complex frequency. | ||
# simplify naming convention for the container for D_k(w) | ||
""" | ||
function tabLeg!_NEUTRAL(omg::Complex{Float64}, struct_tabLeg::LegendreFHT) | ||
# Simplify naming convention for the container for D_k(w) and Q_k(w) | ||
tabDLeg = struct_tabLeg.tabDLeg | ||
tabQLeg = struct_tabLeg.tabQLeg | ||
tabPLeg = struct_tabLeg.tabPLeg | ||
|
||
# keep only the real part of omg | ||
# Keep only the real part of omg | ||
romg = real(omg) | ||
|
||
# initial values of Q_k(w): | ||
|
||
# Initial values of Q_k(w): | ||
if romg == 1.0 | ||
# unique case, requires special seed values | ||
# Unique case, requires special seed values | ||
val_0_Q = 0.0 + 0.0im | ||
val_1_Q = 2.0 + 0.0im | ||
|
||
else | ||
# initial value for k = 0. ATTENTION, not to forget the `abs'. | ||
val_0_Q = log(abs(1.0-romg)) - log(abs(-1.0-romg)) + 0.0*im | ||
|
||
# initial value for k = 1. | ||
val_1_Q = 2.0 + romg*val_0_Q + 0.0*im | ||
# Initial value for k = 0 | ||
val_0_Q = log(abs(1.0 - romg)) - log(abs(-1.0 - romg)) + 0.0im | ||
# Initial value for k = 1 | ||
val_1_Q = 2.0 + romg * val_0_Q + 0.0im | ||
end | ||
|
||
# now the tabQLeg | ||
tabQLeg!(omg,val_0_Q,val_1_Q,tabQLeg) | ||
# Compute tabQLeg | ||
tabQLeg!(omg, val_0_Q, val_1_Q, tabQLeg) | ||
|
||
# compute P_k(w) if needed: | ||
# Computing the value of H[-1,1,Re[w]]. ATTENTION, the argument must be real. | ||
# Compute P_k(w) if needed | ||
val_H = Heaviside(romg) | ||
|
||
# We are outside of the interval [-1,1]: no need to compute P_k(w) | ||
if (val_H == 0.0) | ||
|
||
# compute D_k(w): | ||
|
||
# Loop over the Legendre indices | ||
for k=1:(struct_tabLeg.Ku) | ||
|
||
# Filling in tabDLeg | ||
# We are outside of the interval [-1,1], no need to compute P_k(w) | ||
if val_H == 0.0 | ||
# Compute D_k(w) | ||
for k in 1:struct_tabLeg.Ku | ||
tabDLeg[k] = tabQLeg[k] | ||
end | ||
|
||
# We are within the interval [-1,1]: we need to compute P_k(w) | ||
else | ||
# Initial values of P_k(w) | ||
val_0_P = 1.0 + 0.0im # Initial value for k = 0 | ||
val_1_P = romg + 0.0im # Initial value for k = 1 | ||
|
||
# Initial values of P_k(w): | ||
# Initial value for k = 0. | ||
val_0_P = 1.0 + 0.0*im | ||
|
||
# Initial value for k = 1. | ||
val_1_P = romg + 0.0*im | ||
tabPLeg!(omg, val_0_P, val_1_P, struct_tabLeg.Ku, tabPLeg) # Compute tabPLeg | ||
|
||
tabPLeg!(omg,val_0_P,val_1_P,struct_tabLeg.Ku,tabPLeg) # Computing the tabPLeg | ||
|
||
# Computing D_k(w): | ||
|
||
# Loop over the Legendre indices | ||
for k=1:(struct_tabLeg.Ku) | ||
|
||
# Filling in tabDLeg | ||
tabDLeg[k] = tabQLeg[k] + im*pi*val_H*tabPLeg[k] | ||
# Compute D_k(w) | ||
for k in 1:struct_tabLeg.Ku | ||
tabDLeg[k] = tabQLeg[k] + im * π * val_H * tabPLeg[k] | ||
end | ||
|
||
end | ||
end |
Oops, something went wrong.